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Abstract 

We consider issues associated with the Lagrangian characterisation of flow structures arising 
in aperiodically time-dependent vector fields that are only known on a finite time interval. A 
major motivation for the consideration of this problem arises from the desire to study transport 
and mixing problems in geophysical flows where the flow is obtained from a numerical solution, 
on a finite space-time grid, of an appropriate partial differential equation model for the velocity 
field. Of particular interest is the characterisation, location, and evolution of transport barriers 
in the flow, i.e. material curves and surfaces. We argue that a general theory of Lagrangian 
transport has to account for the effects of transient flow phenomena which are not captured 
by the infinite-time notions of hyperbolicity even for flows defined for all time. Notions of 
finite-time hyperbolic trajectories, their finite time stable and unstable manifolds, as well as 
finite-time Lyapunov exponent (FTLE) fields and associated Lagrangian coherent structures 
have been the main tools for characterizing transport barriers in the time-aperiodic situation. 
In this paper we consider a variety of examples, some with explicit solutions, that illustrate, in 
a concrete manner, the issues and phenomena that arise in the setting of finite-time dynamical 
systems. Of particular significance for geophysical applications is the notion of flow transition 
which occurs when finite-time hyperbolicity is lost, or gained. The phenomena discovered and 
analysed in our examples point the way to a variety of directions for rigorous mathematical 
research in this rapidly developing, and important, new area of dynamical systems theory. 

1 Introduction 

Organised or 'coherent' structures in fluid flows have been a subject of intense study for some 
time, especially since the seminal paper of Brown and Roshko ( 14 ). The dynamical systems 
approach to the Lagrangian aspects of fluid transport, which became widespread in the 1980's 
and 90's, has provided a variety of techniques for determining the existence and quantifying 
'organised structures' in fluid flows. Hyperbolic trajectories and their associated stable and 
unstable manifolds have provided one approach to this problem, in both the periodic and 
aperiodic time dependent settings, that dates back to the beginning of studies of 'chaotic 
advection' in fluid flows (JBH] H \M EH HI])- More recently, the notion of 'Lagrangian 
coherent structure' (henceforth LCS) derived from finite-time Lyapunov exponent (FTLE) fields 
has provided another means of identifying coherent flow structures in fluid flows which can be 
used in Lagrangian transport analysis (42 37, 38, 78, 59 ). The purpose of this paper is to 
compare the methods based on determination of stable and unstable manifolds of hyperbolic 
trajectories with LCS's derived from FTLE's as techniques for uncovering organised structures 
in fluid flows and quantifying their influence on transport. 



We begin in Section [2] by reviewing some theoretical issues associated with Lagrangian 
transport analysis in time-dependent vector fields defined over a finite time interval. We also 
and take the opportunity to clarify a number of misconceptions that have arisen in the literature 
concerning the applicability of hyperbolic trajectories and their stable and unstable manifolds 
in analysing Lagrangian transport in fluid flows, especially with respect to their comparison 
with LCS's. This will naturally lead to the issue of a relationship between the stable and 
unstable manifolds of hyperbolic trajectories and LCS's. 

The purpose of this paper is to compare the methods based on detection of stable and 
unstable manifolds of hyperbolic trajectories with LCS's derived from FTLE fields as techniques 
for uncovering organised structures in unsteady fluid flows. We will particularly focus on the 
performance and applicability of these techniques in flows undergoing transitions associated 
with a loss or gain of finite-time hyperbolicity by some trajectories. An understanding of 
this relationship is essential for understanding the role that each of these structures plays in 
Lagrangian transport. Both methods can have drawbacks as tools for diagnosing the finite- 
time Lagrangian flow structure, In Section [3] we consider a series of examples which aim at 
providing a guide for choosing the most suitable technique for a particular application. We 
begin the discussion by studying a one-dimensional non-autonomous system which can be 
solved analytically and which provides a good illustration of issues concerning the finite-time 
hyperbolic trajectories and FTLE fields in higher dimensions. The subsequent examples of 2D 
non-autonomous systems are chosen to highlight various properties and problems arising in the 
invariant manifolds and FTLE analysis. 

We summarise our findings in Section [4] where we also discuss a number of outstanding 
problems. The Appendices contain a number of technical details and definitions, as well as a 
discussion of some important facts necessary for computation of finite-time stable and unstable 
manifolds. 

2 Some Theoretical Background and Questions 

In this section we describe some of the relevant theoretical issues related to hyperbolic trajec- 
tories and their stable and unstable manifolds and LCS's. This will serve to highlight some 
practical issues arising from applications and computation, as well as the need for further 
theoretical and computational developments. We will not go into great detail in describing 
the theoretical results and computational methods since they are already covered in numerous 
papers in the literature; relevant references will be provided wherever appropriate in the discus- 
sion. Rather, we will discuss ideas and concepts and provide a guide to the existing literature. 
In order to achieve a relative self-containment of the following discussion, we also provide a 
number of important definitions in the Appendix A in order to make this discussion easier to 
follow. 

The notion of hyperbolicity of a trajectory has been around for some time. It is particularly 
worth remembering in the context of the present discussion that hyperbolicity is not dependent 
on the nature of the considered time dependence (although continuity in time, which is also our 
operating assumption here, eliminates many technical issues). In particular, if hyperbolicity is 
determined by Lyapunov exponents ([51]) or exponential dichotomies ([19]), then the nature 
of the time dependence, e.g. periodicity, quasiperiodicity, or aperiodicity plays no role in any 
of these definitions (and equivalence between these definitions is considered in [24]). Once a 
hyperbolic trajectory is located, then the stable and unstable manifold theorem for hyperbolic 
trajectories immediately applies, and this is also independent of the nature of the time depen- 
dence. It can be verified that the statement of this theorem is also independent of the nature 
of the time dependence by examining, for example, its proof in the classic ordinary differential 
equations textbook of Coddington and Levinson [TH]. Additional resources on the stable and 
unstable manifold theorem for arbitrary time dependence can be found in |2lT 1461 151] . 

Of course, a central issue in practical applications is the location of hyperbolic trajectories 
in aperiodically time dependent velocity fields. Historically, there have been many algorithms 
for finding equilibrium points (stagnation points) of steady velocity fields and periodic orbits of 
time-periodic velocity fields, but relatively little work had been done on algorithms for finding 
hyperbolic trajectories of aperiodically time dependent velocity fields (and quite a few new is- 
sues arise, in comparison to the issues associated with steady and time periodic velocity fields, 
which we will mention below). An algorithm for determining hyperbolic trajectories in arbi- 
trary unsteady flows was given in [^5] and further refined in [49 1 161] . This technique is based 
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on an iterative method denned on a space of 'paths' and, provided it converges, is guaranteed 
to yield a hyperbolic trajectory on a specified time interval which is bounded in most practical 
applications. (The 'finiteness' of the considered time interval brings up yet another technical 
issue that we will shortly address.) The iterative algorithm requires an initial 'guess' in the form 
of a C 1 path defined on the appropriate time interval. It is important to stress here that such 
a path need not be a trajectory of the velocity field. We provide a few more details regarding 



some necessary properties of the initial guess in the Appendix A (cf Definition A. 5 and remarks 



after Definition A. 6 1. The initial guess is often chosen to be a path of hyperbolic instantaneous 
stagnation points (ISPs, cf ( |83| ), Appendix A). This particular choice of the initial path has lead 
to numerous misleading and incorrect statements in the LCS literature related to the notion of 
"Galilean invariance" and the nature of this algorithm ( [581 1591 178] ). Galilean transformations 
consist of spatial translations, time translations, shear transformations, reflections, and rota- 
tions. Paths of ISPs are not, in general, particle trajectories. This has been a known fact in 
the fluid dynamics community for some time, and a simple proof can be found, for example, an 
appendix in [35]. Clearly, ISPs are not invariant under Galilean transformations. However, it is 
well-known in the dynamical systems community that trajectories are invariant under Galilean 
transformations (i.e. a trajectory maps to a trajectory under a Galilean transformation) and 
hyperbolic trajectories to which the iterative algorithm converges are likewise invariant under 
Galilean transformation^] Consequently, the fact that a non-Galilean invariant path is used as 
an initial guess for the iterative algorithm is irrelevant since if the algorithm converges, it yields 
a hyperbolic trajectory, which is manifestly Galilean invariant. Likewise, since the stable and 
unstable manifolds of a hyperbolic trajectory are, by definition, composed of trajectories, they 
are also Galilean invariant. The importance of Galilean invariance to specific oceanographic 
investigations is another matter entirely. Oceanographers require a fixed reference frame to de- 
scribe the ocean through measurements and grid based computations. In the chosen frame, the 
behavior and stability of ISPs have historically played an important role in describing observed 
Eulerian flow structures. While ISPs may bear little relation to particle trajectories, we believe 
that dismissal of their utility on the grounds of not being Galilean-invariant is unjustified. 

Now we return to a more serious issue. Hyperbolicity, and therefore hyperbolic trajectories 
and their stable and unstable manifolds are 'infinite-time objects'. More precisely, hyperbol- 
icity of a trajectory is determined on the basis of the asymptotic behaviour of neighbouring 
trajectories in the infinite time limit. The stable and unstable manifolds associated with a hy- 
perbolic trajectory are proven to exist via a fixed point, or iterative, argument where the limit 
as time goes to either positive or negative infinity is taken. If the velocity field is aperiodic in 
time, and it is obtained from the output of a numerical computation, then we have knowledge 
of the velocity field only on a finite time interval. This fact creates a host of new problems in 
applying the 'traditional' dynamical systems approach to fluid transport. The main difficulty 
in the 'finite-time' description of Lagrangian transport stems from the fact that the dynamical 
systems theory is generally concerned with the 'long time behavior' of systems of ODE's (many 
of these problems are discussed in [86] and [62]). In particular, the standard definitions of 
hyperbolicity of trajectories do not apply to velocity fields that are only known on a finite time 
interval (henceforth finite-time velocity fields). 

The subject of 'finite-time dynamical systems theory' gives rise to many new issues that 
require new theoretical and computational results. These are discussed in [861 162] . There have 
also been a number of mathematical papers developing various aspects of this subject in recent 
years (|28l \8\). The 'finite-time' framework is intrinsically dependent on the time interval one 
considers in the analysis and the implications of non-uniqueness associated with this setting 
has been discussed in numerous papers, see, e.g., [B(j] E3 I3S1 ESI ED- ^ n particular, in the 
context of finite-time dynamical systems, hyperbolicity of a trajectory is defined over a finite 



time interval (cf Definitions A. 4 and A. 11 in the Appendix A) and the stable and unstable 
manifolds associated with the trajectory no longer have a lower dimension than the underlying 
phase space (cf Appendix [B] and [28]). Consequently, a trajectory which is hyperbolic over 
some time interval (in the finite time sense) may not be hyperbolic over a longer time interval. 
In other words, given that a < b < c < d, it is possible for a trajectory to possess finite- 
time hyperbolic characteristics on all intervals contained in I a b = [a, 6], and then lose such 
characteristics on some intervals contained in he, possibly regaining the finite-time hyperbolic 
properties for all intervals contained in Id- We refer to such a scenario as a 'loss' and a 



1 The Galilean invariance of hyperbolic trajectories is proven in [45| for hyperbolicity determined with exponential 
dichotomies. 
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subsequent 'gain' of finite-time hyperbolicity and point out that one cannot pin these transitions 
to a particular time instant. Purists in dynamical systems theory may immediately object by 
saying that hyperbolicity is a notion that only has meaning for trajectories defined for all 
time. According to the traditional definition, this is certainly correct. However, applications 
to transport in velocity fields defined for finite time have motivated this new definition of 
hyperbolic-like properties over a finite-time (i.e. the finite-time hyperbolicity) and the notion of 
loss or gain of (finite-time) hyperbolicity has proven useful for describing the transient behavior 
of a number of time dependent structures in oceanographic flows. We will discuss examples of 
simple flows whose transitions are induced by the loss (or gain) of finite-time hyperbolicity in 
pT2l pT6] and pTT) 

In any case, it is important to realise that all of the finite-time dynamical systems notions 
that we mentioned above are trajectory based. That is, the finite-time hyperbolic trajectories 
are indeed trajectories and material curves contained in their finite-time stable and unstable 
manifolds are barriers to transport (see also Appendix Their usefulness for applications 
derives solely from their ability to explain new phenomena in applications, and this is assessed 
in the context of specific applications. 



We now turn to another technique used in the finite-time transport analysis which is based 
on determination of the so-called Lagrangian coherent structures (LCS) from finite-time Lya- 
punov exponent fields (FTLE). Lyapunov exponents are quantities associated with trajectories 
that are obtained as infinite time limits. For an n-dimensional continuous time dynamical sys- 
tem a trajectory has n Lyapunov exponents - one associated with a direction tangent to the 
trajectory (which is always zero) and n — 1 Lyapunov exponents associated with the remaining 
directions. The Lyapunov exponents are measures of the growth of infinitesimal perturbations 
in these directions, i.e. growth rates of the linearized dynamics about the trajectory (cf Ap- 
pendix DQ. Of particular interest is the maximum Lyapunov exponent since the existence of a 
single positive Lyapunov exponent indicates that the trajectory is unstable. The fundamental 
theorem on the existence of Lyapunov exponents is expressed by the Oseledec multiplicative er- 
godic theorem ([67])- There are many excellent references on Lyapunov exponents that describe 
their properties ( \51\ 1561 157p and algorithms for their computation ([23, 22, 24, 35, 32 ). 

In the infinite-time setting, Lyapunov exponents are one measure of the hyperbolicity of 
a trajectory. If a trajectory has nonzero Lyapunov exponents (with the exception of the zero 
exponent associated with the direction tangent to the trajectory), it is said to be hyperbolic 
f |51|). Finite time Lyapunov exponents are obtained by computing the same quantities, but 
restricting the computation to a finite time interval, rather than taking the limit as the time 
goes to positive infinity (for forward time Lyapunov exponents) or minus infinity (for backward 
time Lyapunov exponents^ Clearly, one would like to know the length of the time interval 
on which they must be computed so that they are "close" to the infinite time limit. Some 
interesting arguments are given in [341 129] which indicate that the rate of convergence may 
be quite slow. The FLTE technique]^] is not immune to the non-uniqueness issues arising in 
the finite time setting mentioned earlier. These are highlighted by the fact that for any time 
instant in the considered time interval I one can compute a whole family of FTLE fields. We 
discuss implications of this fact in the following sections. 

For each time instant t n within the considered (or available) time interval I, forward FTLE 
fields are obtained by computing the forward Lyapunov exponents of the trajectory starting at 
that initial condition at t n in a chosen grid for the length T of time available (and computable) 
and colour coding the initial condition according the the magnitude of the largest FTLE (e.g. 
bright colors for large values, light colors for small values). By performing such a computation 
for an ordered sequence of 'observation times', {t n } ne z,t n £ I, one can examine the spatial 
evolution of the structures exhibited by the forward FTLE fields in time. Clearly, backward 
FTLE fields can also be computed by reversing the direction of time. Note here that for any 
t n in such a sequence it is possible to compute an FTLE field for any T such that t n + T £ I . 
It is often not obvious which length of the integration time interval T should be chosen in such 
computations especially when the structure of the resulting FTLE fields varies significantly for 



2 We note that in the literature the notion of a "direct Lyapunov exponent" (DLE) has been introduced (|37| h 
This has created some confusion in the literature in the sense that the acronyms "FTLE" and "DLE" are used 
somewhat synonymously. In recent years the consensus has become that there is no substantive difference between 
the two notions and "FTLE" has now returned to being the accepted acronym (e.g., see 78, 76 77 59 ). 

3 We note that in much of the literature concerning FTLEs, the phrase refers to the maximum FTLE. 
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different values of T. We discuss these issues in most of the examples presented in ^3] 

Since Lyapunov exponents are a measure of the (linearized) growth rates of a set of or- 
thogonal directions perpendicular to the tangent vector to a trajectory, FTLE fields have been 
more physically referred to as "stretching fields']^] Numerous groups have computed FTLE 
fields over the years in the context of fluid transport (e.g., |69l 1701 183] ) and have noted that 
these fields appear to exhibit a great deal of structure. A more precise quantification of such 
structures have led to the notion of LCS ([36l [42l [3Tj [38] [39l [T8] [59]). In particular, since 
FTLE's are a measure of separation of nearby trajectories after some finite-time, regions of 
high values for the maximal FTLE would seem to be likely candidates for regions containing 
hyperbolic trajectories and their stable and unstable manifolds. Heuristic arguments support- 
ing this assertion are given in the aforementioned references, and will not be reproduced here. 
Rather, in this paper we will focus upon the assumption that "maxima" of the FTLE fields 
are "approximations" to the unstable manifolds of hyperbolic trajectories (forward time FTLE 
fields) and unstable manifolds of hyperbolic trajectories (backward time FTLE fields). We have 
put the word maxima is quotes since this notion needs careful consideration. This was done in 
|78| via the notion of a ridge curve of an FTLE field. Roughly speaking, a ridge curve has the 
property that moving transverse to the direction tangent to the curve corresponds to moving 
to a lower value of the FTLE. Precise definitions are given in 78 where ridges of the FTLE 
field are taken as the definition of LCS. This raises the question of precisely how "Lagrangian" 
are LCS's? In general, they are not material curves, and therefore not necessarily barriers to 
transport. In the following sections we will demonstrate this with several examples designed 
to highlight different aspects of the problem. Nevertheless, certain segments of an LCS may 
be "close" to a barrier to transport in the sense that the flux across the curve may be small. 
This issue was carefully considered in [7S]. However, the extent to which LCS's are barriers 
to transport must be assessed after they are computed. The stable and unstable manifolds of 
finite time hyperbolic trajectories are a priori barriers to transport since they are computed as 
curves of fluid particle trajectories. 

We remark that a possible misconception that has appeared in several places in the LCS 
literature is that the concept of invariant manifold is somehow either not well defined or applica- 
ble or easily interpretable for time-dependent flows, where the time dependence is not periodic 
f [421 1371 1381 1761 '58 ). In particular, this point has been emphasized in the finite time dynam- 
ical systems context. While the approach to Lagrangian transport based on finite-time stable 
and unstable manifolds of finite-time hyperbolic trajectories certainly requires more complex 
algorithms and computational techniques, the results, being trajectory based, are certainly un- 
ambiguous (in that sense) and the value of the approach can only be assessed in its ability to 
explain Lagrangian transport phenomena. Towards this end we note that [BIJ] utilises a finite 
time, realistic velocity field obtained from a data assimilating oceanographic model (DieCAST) 
that uses finite time hyperbolic trajectories and their (non-unique) stable and unstable mani- 
folds to give the first Lagrangian characterization of a salinity front in the Mediterranean Sea 
and provide and explanation and characterization of the notion of 'leakiness' of the front. Of 
course, the finite time issues mentioned above do require careful consideration in the context 
of specific applications. It is incorrect to think that the LCS approach has somehow "solved" 
this problem. 

A broader issue here, which keeps recurring throughout the following discussion, concerns 
the problem of description of the Lagrangian structure of a time-dependent flow in a way which 
would allow for a meaningful finite-time Lagrangian transport analysis. It is well known that 
in order to establish the existence of, for example, a transport barrier (i.e. a flow-invariant, 
Lagrangian structure) in the non-autonomous case, one requires non-local (in time and space) 
information about the governing flow. As already pointed out, the finite-time notions discussed 
above may provide ambiguous diagnostics due to their potential sensitivity to the time-interval 
chosen for extracting the relevant information. Consequently, it seems crucial for the devel- 

4 As we have noted, FTLE's are a measure of the growth of "infinitesimal perturbations" to a given trajectory, i.e. 
growth rates of the linearized dynamics about a trajectory. Finite size (or "scale") Lyapunov exponents (FSLE's) are 
a technique to analyze the growth of "finite perturbations" to a given trajectory. Alternatively, FSLE quantify the 
relative dispersion of two particles, as discussed in [9]. In 9, 53, 48 25 31 26 Lagrangian structures are identified 
using FSLE's. The maxima of the FSLE fields look very much like the maxima of FTLE fields and bear a striking 
resemblance to the stable and unstable manifolds of hyperbolic trajectories. However, it must be emphasized that 
FSLE's are a non-rigorous numerical technique and, despite the strong numerical evidence, there are no theorems 
relate the results of the calculations to Lagrangian transport barriers. Much like the case with FTLE's, this must be 
assessed "after the fact" . 
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opment of a general theory of finite-time transport in aperiodically time-dependent velocity 
fields to understand and properly describe transient flow phenomena. Undoubtedly, this task 
requires development of tools which would adequately capture the finite-time flow properties. 
The examples discussed in the next section highlight a number of important points regarding 
the techniques of invariant manifolds and FTLE fields: 

(1) One can obtain a good agreement between the ridges of the FTLE fields (i.e. the LCS) 
and the finite-time stable/unstable manifolds of distinguished hyperbolic trajectories in 
sufficiently 'well-behaved' flows, 

(2) Both approaches may provide non-unique results, particularly in flows undergoing tran- 
sitions (discussed later), and their interpretation may require a subjective interpretation. 
The main drawback affecting the invariant manifold computations lies in identifying the 
appropriate hyperbolic trajectory (i.e. the Distinguished Hyperbolic Trajectory) used for 
'seeding' the finite-time stable and unstable manifolds. The main drawback affecting the 
FTLE technique stems from the fact that it is a function of trajectory separation which 
depends, in general, on the time interval chosen for assessment of such a measure. Conse- 
quently, in flows undergoing transitions it is often difficult to decide which time interval is 
most suitable for assessing the (non-local) flow structure. Moreover, there is no guarantee 
that the time evolution of the ridges of locally strongest separation is continuous in time. 



3 Tests 

In this section we analyse a wide range of example flows for which both the FTLE fields and 
the appropriate invariant manifolds are computed. We then analyse and compare the infor- 
mation about the Lagrangian flow structure obtained from computing the backward/forward 
FTLE maps, and the information obtained from computing the unstable and stable manifolds 



of certain Distinguished Hyperbolic Trajectories (cf Definition A. 6 1 in these flows. The al- 
gorithms used for computing the DHTs and their manifolds, based on the ideas described in 
[HO BH |H31 ED , were developed in MATLAB. The FTLE computations are performed also in 
MATLAB using an implementation of methods described in |37l 1781 l76l 177] . We also compare 
our results with the LCS MATLAB Kit v. 2. 3, developed in the Biological Propulsion Labora- 
tory at Caltech, which is available online [201 ■ In the case of the LCS MATLAB Kit, several 
minor modifications were introduced in the code in order to enable FTLE computations from 
analytically defined vector fields. 

All the examples considered here are based on analytically defined velocity fields. While 
the resulting flows are certainly not sufficiently complex to be of importance in practical ap- 
plications, they provide an easily reproducible testbed for our analysis. 



3.1 ID non-autonomous configuration 

We consider first a one-dimensional, non-autonomous ODE which can be solved analytically, 
and which illustrates in the simplest possible setting a number of issues which are important 
in the following sections. Based on three related examples, we highlight potential difficulties 
when trying the uncover the structure of a non-autonomous flow using the finite-time Lyapunov 
exponents, or when trying to identify some 'special' trajectories which play an important role 
in organising the global dynamics. Of course, in such a setting there are no non-trivial invariant 
manifolds in the (non-autonomous) flow. However, one can consider the ID geometry discussed 
below to represent some aspects of transverse dynamics in the neighbourhood of an invariant 
manifold in a higher-dimensional flow; in fact, we use this analogy in §3.2.7| Here, we are 
particularly interested in the properties of the FTLE maps and their relationship to during 
certain flow transitions characterised by changes of finite-time stability properties of some 
distinguished trajectories in the flow. 



Consider a one-dimensional, non-autonomous dynamical system given by 

x = x(a(t) - x 2 ), x,teR, (1) 

where a(t) is a prescribed function of time. In the autonomous configuration, with a — const. < 0, 
the trivial solution x = 0, representing the only fixed point in the flow, attracts all trajectories 
as t ^ oo. When a = const. > 0, there are three fixed points in the flow: x\ — 0, and 



G 



Scenario I Scenario II Scenario III 




Fi gUT6 1; Geometry of the one-dimensional flows with the time dependence induced by cr(t) characteristic of the 
three scenarios considered in |3.1| The trajectories 71 (t), 72 (t) are distinguished in the sense described in appropriate 
sections. Analysis of these flow structures using the FTLE technique are summarised in figures [2] [3] and [4] 



X2,3 = ±yo\ It can be easily checked by examining the linearisation of |TJ about these points 
that xi is an unstable hyperbolic fixed point and £2,3 are stable hyperbolic fixed points. 

When da/dt 7^ 0, it is more convenient to consider the resulting dynamics in the extended 
phase space, spanned by {e x ,e t }, with coordinates (x,t). We note here that jlj) is, in fact, a 
Bernoulli equation with solutions given by the family 



x{t;x ,t ) 2 = — . (2) 



X J t 



It can be easily verified using |2| that x(to, xo,to) ~ xo- For any trajectory x(t,xo,to), given 
by (pj, we can consider a perturbation, x(t, xo + So, to), with So <C 1, so that the growth of the 
perturbation after time T is given by 



S(T, 5o,x , t ) = \x(t + T, x ,t ) — x(t + T, x + S , t )\ 



dx(t + T, s,t ) , s - v . _> 



6o + 0(V) 
(3) 

Thus, since the solutions Q are continuous, the growth of an infinitesimal perturbation intro- 
duced at (xo,to) after time T is given by 



A(T,x ,t ) = lim 



5(T, S , t ) 



-2 f* 0+T a(s)ds 



^o^o So 



-2/*° + T a( S )ds , „ 2 



3/2 



(4) 



We note further that Q is related to the ID finite time Lyapunov exponent \t{xo^o) at time 
to via 

\ T {xo,to) = lnA(T,xo,to), (5) 

which is computed over the time interval T, (see the Appendix for a more general formulation). 

Note that even if solutions satisfying a given system are only known numerically, an estimate 
on the separation rate of trajectories which were initially infinitesimally close can be obtained 
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b) c) 




Figure 2: (a-d) ID FTLE fields, \ T {x, to), for the flow |Tj with a{t) given by | |l6| which is characteristic of Scenario 
I discussed in j |3.1| The finite time Lyapunov exponents, Ay, are computed over different time intervals of length T. 
In this configuration, there are three 'distinguished' trajectories in the flow, 71,2(4) (cf JlOjl) and x = 0, which play 
an important role in organising the dynamics (blue curves; left column), (b-c) Backward FTLE field computed, using 
J4| and ^5), at (b) t = 5 and (c) t = — 1 with different values of the integration parameter T. Note that the maxima 
of the FTLE fields (i.e. the LCS) vary with T, and that they do not coincide with the location of 71,2(4 = — 1) in the 
transition phase (e.g. (c)), regardless of the value of T. See text for a discussion, (d-e) Forward FTLE field computed 
for the same flow at (d) t = and (e) t = — 10 with different values of the parameter T. 



via finite differences. Therefore, At can be estimated for any flow defined by sufficiently smooth 
velocity field on some time interval /. Consequently, the map 

R3 x -> At(Mo) 6 JR, to+TelCM, (6) 

can be used, in principle, as a straightforward diagnostic tool for uncovering time-dependent 
flow structures characterised by locally strongest separation of nearby trajectories. Note how- 
ever, that at any time to during the flow evolution one can construct the whole family of 
FTLE fields {Xt(x, to)}T+t ei which generally results in a non-uniqueness of the computed 
diagnostic. The ambiguities associated with choosing the 'right' FTLE map from the family 
{Xt(x, to)}T+t £i which 'best' describes the flow structure at a given time are especially evi- 
dent in analysis of flows displaying transient phenomena. We recall that this problem is not 
restricted to the FTLE method. In particular the techniques, mentioned in ^2] based on iden- 
tification of the so called 'distinguished hyperbolic trajectories' and their invariant stable and 
unstable manifolds suffer from similar limitations in the case of flows defined on a finite time 
interval. We analyse these issues further below based on three different scenarios of evolution 
of the one-dimensional flow Q, characterised by different types of time dependence induced 
by the form of cr(t). Clearly, the dimensionality of the problem does not allow for existence 
of any non-trivial invariant manifold of a hyperbolic trajectory. Nevertheless, the discussed 
examples serve to highlight some important consequences of flow transitions (specified below) 
on the computed FTLE fields and their relationship to some (possibly non-unique) 'special' 
trajectories in the space of solutions of Moreover, we will show that the non-uniqueness 
of the FTLE diagnostic may lead to detection of 'ghosts' or 'premonitions' of flow structures 
associated with the future, or past, stability properties of such 'special' trajectories. We will 
later return to these examples in |3.2.2| in the context of locally transverse dynamics in a neigh- 
bourhood of a stable or unstable manifold of a hyperbolic trajectory in the 2D non-autonomous 
case. 
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Scenario I: < a(t) < oo. 

With the above constraints imposed on cr(t), the trivial solution, x(t) = 0, of |l]) is (finite- 
time) unstable on any time interval I = [t a ,if>] G H in the sense that for each nonempty, 
bounded set xi 9 there exists a trajectory, x(t, xo, to), with xo G Xi, to G I, such that 

^|a:(*,x ,to)| > 0, Vt 6 J. (7) 

A more general definition of instability of a trajectory in a non-autonomous dynamical system, 
which we do not require here, can be found, for example, in [55]. It can be easily verified that 
|7f is satisfied on x{t) =0 over any time interval / C IR by noticing that 



d_ 
dt 



\x(t,x ,t ) 2 J *0 Jt 
which implies that |7| is satisfied at least for 



e- 2f ^^ de dk + x 2 ), (8) 



2 "mint- - . 

x ° 1 _ e -2 CTmax (t b -t a ) ( - e 2 CTmin (t b -i a ) _ x y W 

We note further that there are two 'distinguished' trajectories in the space of solutions of 
|TJ given by 

71.2W 2 = — , (10) 



which have the property that any trajectory of |T]) x(t, xo,to), xq > is 'attracted' (in the sense 
we specify below) towards 71 (t) and any trajectory x(t, Xo, to), xo < is 'attracted' towards 
72 (t). There are two different notions of attraction which we can utilise here. If we rewrite |2| 
as 

2 1 
x(t;x ,to) =- w , (11) 

1 e ~ 2 St a -W d » + _1 2I e- 2 ti CT(a)ds 



< ' 7(*) 

it can be seen that the following are true (when < a(t) < oo) 



t — ►oo 



lim [x(t, xo, to)- 71 (*) =0, Vx o <0,t elR, (12) 



Km \x(t,xo,to)-j2(t) =0, Vxo>Q,t €]R, (13) 



and 



lim ( x(t,x ,t ) 2 - 7(t) 2 ) = 0, Vi <0,fefi, (14) 

to— ►— oo \ / 

lim (a;(t,a;o,to) 2 -7(t) 2 ) = 0, Va;o<0,tG7R. (15) 
to— ►— oo \_ y 

Since we intend to minimise the amount of mathematical formalism here, we just remark 
that the property ( |12[) implies that 71 (t) is forwards attracting (and Lyapunov stable) within 
xo G (—oo,0) and ( 14 1 implies that it is pullback attracting within xo G (— 00, 0). Similarly, 
72(t) is both forwards and pullback stable within xo G (0, 00). A more formal introduction to 
the stability and bifurcation phenomena in non-autonomous dynamical systems can be found 
in [551 [54 , 52 , 28 72 73 . Pullback convergence is useful in constructing limiting sets, such as 
the distinguished trajectories in our ID toy example, provided that the flow is defined on the 
negative half-line (— 00, t*], t* > —00. Otherwise, we cannot uniquely define a distinguished 
trajectory. We will see in the next example that these two notions are not necessarily equivalent 
in the non-autonomous case. 

We can now examine the one-dimensional FTLE fields, At(:e, to), associated with scenario 
I which are obtained from (5| and |4| for different lengths of the integration time interval, T. 
The results shown in figure 2] were computed for a sigmoidal function 

ff(t) = -(atan(10(t + 4)) + 7r/2 + 0.01), (16) 

7T 
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b) c) 




Figure 3: (a-d) ID FTLE fields, X T (x,t ), for the flow |T} with a(t) given by \20\ which is characteristic of 
Scenario II discussed in |3.1| the fields, Ay, arc computed over different time intervals of length T. (b-c) Backward 
FTLE field computed, using fil and JSl, at (b) t = 5 and (c) t = —8 with different values of the integration 
parameter T. In this configuration there distinguished trajectories 71,2 W (cf l |10| ) dominate the flow structure after 
the transition when the trivial solution becomes unstable. Note that for sufficiently large values of the integration 
parameter T the maxima of the FTLE fields detect 'ghosts' of the past stability of the trivial solution and not 
the situation at the time of computation t. See text for a discussion, (d-e) Forward FTLE field computed for the 
same flow at (d) t = — 12 and (e) t = with different values of the parameter T. The trivial solution x = is 
globally attracting in the sense of (|17|l on any time interval / = (— oo,t^_],t*_ < t* where t* ~ —4.105. Note that, 
when computed over sufficiently long time intervals, the FTLE fields detect 'premonitions' of the future (finite-time) 
stability properties of the trivial solution (cf (d)) which is repelling (in this case) on any time interval contained in 
/ = (-4.105,oo). 



so that the flow |T| is asymptotically autonomous. 

The top-row insets of figure [2] focus on detection of attracting structures in the (extended) 
phase space of the flow |TJ. Since such structures should be characterised by separation of 
trajectories in backward time, we compute a number of the backward FTLE fields at two 
different times t — 5 (b) and t — — 1 (c). The geometry of the two attracting distinguished 
trajectories 71,2 (i) is marked by the blue curves. Note that the maxima of the FTLE fields (i.e. 
the LCS) vary with T, and that they do not coincide with the location of 71,2(2 = —1) in the 
transition phase (e.g. (c)), regardless of the value of T. The maxima of the forward FTLE fields, 
computed for the same flow at (d) t = and (e) t = —10, are all located at the trivial solution 
x = which is unstable. However, during the flow phase when the unstable trivial solution 
is 'sandwiched' between the two attracting 'distinguished' solutions 71,2, the FTLE field has 
to be computed over sufficiently long time intervals in order to reveal a positive maximum 
(i.e. exponential growth of the infinitesimal perturbation to x = over the considered time 
interval) . 

Scenario (II): lim a(t) < 0, a{t*) = 0, and do /At > 0. 

t — >— 00 

In this situation the trivial solution of jl]), x(t) — 0, is stable (in the pullback sense) on any 
time interval / = (— oo,t!_],tl < t* , i.e. 

lim x(t, x ,t ) = 0, V tel, (17) 

to — > — 00 

and unstable, in the sense Q, on any time interval contained in I = (t*,oo). Note that the 
trajectories 71,2(4) (10 1, which are still solutions of |TJ, are now only asymptotically attracting, 
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i.e. 



i — ►oo 



lim x(t,a;o,to)-7i(*) =0, V x < 0,t G R, (18) 



lim 



(^x(t,x ,to)-j2(t)\ =0, V xo > 0, to G -K, (19) 



but they are not asymptotically pullback attracting. We will loosely refer to t* as the transition 
time, since it corresponds to the boundary of the pullback stability of the trivial solution. 
In figure [3] we analyse the phase-space geometry of the flow ([I]) with a(t) given by 

a{t) = tt/2 + 0.8 ( atan ( 10 (* + 4 )) + °' 8 ) ' ( 2 °) 

which satisfies the constraints characteristic of this scenario and changes sign at t* w —4.105. 
Moreover, such a choice introduces an additional simplification to the problem, making it 
asymptotically autonomous. This configuration makes it easier to observe the emergence of 
an 'attracting' structure developing around the trajectories 71,2 (i) after the transition (see 
figure I3J. The FTLE fields, At, shown in figure [3|b-e) are computed using Q and ([5| at four 
different times and over different time intervals of length T. The examples of the backward 
FTLE fields, computed at (b) t = 5 and (c) t = —8 highlight some typical characteristics of 
this technique when applied to flows with transient phenomena. When computed at times 
after the transition (as in (b)) over sufficiently short time interval lengths T, the maxima of the 
FTLE fields coincide well with the location of the distinguished trajectories (dashed blue lines 
in figure |3]^b) ) . Note, however, that for sufficiently large values of T the maxima of the FTLE 
fields detect 'ghosts' (red) of the past stability of the trivial solution and not the situation at 
the time of computation t. It is worth remembering here that while the geometry of the flow 
trajectories and the transition time is known in the considered example, it may not be at all 
obvious what length of the time interval one should choose when computing FTLE fields for a 
realistic, higher-dimensional geophysical flow. A similar problem might occur when trying to 
identify structures characterised by trajectory separation in forward time via the computation of 
forward FTLE fields. We show examples of such computations for the same flow in figure[3jd,e) 
which are computed at (d) t = —12 and (e) t = with different values of the parameter T. 
As already mentioned above, the trivial solution x — is asymptotically pullback attracting at 
any t contained in I = (— oo,£!l], t*_ < t* ~ —4.105. Therefore, no trajectory separates, in the 
sense Q, from the trivial solution on I. The FTLE fields computed in figure^d) correspond to 
such a situation. However, if one computes the forward FTLE fields at t — — 12 for sufficiently 
large T a sharp positive maximum appears which might be interpreted as a 'premonition' of 
the future (finite-time) stability properties of the trivial solution after the transition. 

Scenario (III): a(t) > for t G [**,***], and a(t) < V t G (-co,**] U [t**,oo]. 

In this configuration the trivial solution is the only 'distinguished' one. It is globally asymp- 
totically pullback stable on any time interval I = (—00, t*],t* < t*, i.e. 

lim x(t,x ,t ) = 0, V teI,x eR, (21) 

tQ — * — 00 

and is globally asymptotically stable on any time interval I = oo),t+* > t** , i.e. 

lim x(t, x , t ) = 0, V t e I,x £ R. (22) 

t — >oo 

However, it can be easily verified by examining (J2J that x(t) = is unstable, in the sense of 
condition on any time interval contained in I = [t*, t**]. 

In order to illustrate the typical properties of the FTLE field in such a case we choose the 
time dependence in the following form 



a(t) = 2(V t2/16 -0.4), 



(23) 



so that t* ~ —3.83 and t** ~ 3.83. In figure [4] we examine the backward (b,c) and forward 
(d,e) FTLE fields for this flow configuration, which are computed for different lengths, T, of 
the time test interval. The trivial solution is unstable on any time interval contained (in this 
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Figure 4: (b-d) ID FTLE fields, X T (x,to), for the flow |T} with a(t) given by ( |2->[ i which is characteristic of 
Scenario III discussed in j |3.1| The trivial solution, x = 0, is asymptotically attracting on the time interval I = R and 
globally pullback stable (see | |21[ l) on any time interval / = (— oo,t^_],t* < t* (in this case t* ~ —3.83; see text). 
The trivial solution is unstable on any time interval contained (in this case) within / = [—3.83, 3.83]. (b-c) Backward 
FTLE field computed, using |4j and J5J|, at (b) t = 3 and (c) t = 8 with different values of the integration parameter 
T. Note that the maxima of the FTLE fields (i.e. the LCS) vary with T and, for sufficiently large T, the FTLE 
fields detect a 'ghost' of the past attracting phase of the trivial solution x = (red curves in (b-c)). See text for a 
discussion, (d-e) Forward FTLE field computed for the same flow at (d) t = —8 and (e) t = with different values of 
the parameter T. Note that at t = —8 (d), when x = is attracting and globally pullback attracting, the FTLE field 
computed over sufficiently long interval T detects a 'premonition' of the future unstable phase of the trivial solution. 



case) within / = [—3.83,3.83]. The backward FTLE fields, \T(x,to), computed at t = 3 show 
a similar behaviour as in figure ^b) except that the magnitude of 'ghost' maximum (red), 
indicating the past attracting properties of the trivial solution, is similar to those computed 
for T — —5 and T = — 10. This simple example indicates the possible problems with inter- 
pretation of the families of FTLE fields at time t, {At (x, t)}T+tei, and the right choice of the 
time integration interval best describing the flow structure at the given time t. The forward 
FTLE computations reveal similar ambiguities when trying to detect structures characterised 
by separating trajectories in forward time. The FTLE field computed at t — —8 (d) with T = 2 
indicates correctly the lack of trajectory separation points. The profile of Xio(x,t = —8) is, 
however, rather broad and one might be tempted to increase the integration time interval T 
in order to obtain a more localised profile. If one then computes the forward FTLE field at 
t = —8 with T — 10, the Xi0(x,t = —8) reveals a positive maximum at x = (red curve in 
(d)) which indicates that the perturbations of the trivial solution will eventually separate with 
a positive At- It is important to understand here that this is not an erroneous result. Indeed, 
we know that the trivial solution is unstable on the time interval / = [—3.83, 3.83] and if one 
follows trajectories from t — —8 to a time contained within this interval this is certainly what 
is going to happen. Moreover, if we follow such trajectories to times beyond I, the positive 
maximum disappears again (e.g. \2o(x,t — —8) in figure|4ji). An important question arises in 
connection to this: Which FTLE fields from the T-parametrised family {Xt{x, t)}T+tei, best 
describes the flow structure at t and how do we recognise that it is not always the field with 
sharpest maxima? 

3.2 Two-dimensional, time-dependent flows 

In the reminder of this paper we consider 2D flows which are defined analytically so that there 
is no additional ambiguity with data handling. In each case we determine the stable and 
unstable manifolds of Distinguished Hyperbolic Trajectories and compare the results with the 
LCS identified from the FTLE maps. 
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3.2.1 Two examples of dynamical systems where the Lyapunov exponents 
of every trajectory are equal 

In this section we point out two situations where the Lyapunov exponents of every trajectory 
are equal. Interestingly, the two flows are, in some sense, almost exact opposites in terms of the 
complexity of the dynamics that they exhibit. The first example is the velocity field due to a 
linear, time-dependent straining flow defined on the plane. In this case we can compute the form 
of the Lyapunov exponents explicitly, and as a result it is evident that the Lyapunov exponents 
do not depend of the initial condition of the trajectory, which implies that they are identical for 
all trajectories. In this case the FTLE field reveals no LCS's, for any time over which the FTLE 
field is defined. The second example is the Arnold cat map, defined on the torus (i.e. doubly- 
periodic boundary conditions). It is a linear map, defined on a nonlinear, bounded phase space 
(i.e. the torus). The Lyapunov exponents for every trajectory can be computed explicitly, and 
linearity of the map implies that all exponents are equal. Hence, also in the Arnold cat map 
case the FTLE fields reveal no LCS's. Contrasting these two examples is interesting. Neither 
example has LCS's as diagnosed by the FTLE field (although the phase space of each does 
have hyperbolic trajectories with stable and unstable manifolds), and the velocity field given 
by the linear, time-dependent straining flow has "simple" trajectories, while the trajectories 
exhibited by the Arnold cat map are "extremely" chaotic. We will now describe each of these 
examples in more detail, and in the process provide more background and justification for these 
statements. 



Linear, time-dependent strain: We consider here the simplest class of incompressible 
2D flows, defined for all t £ M, which possess a DHT in the sense of [45]. The flows are trivial, 
time-dependent extensions of the linear steady strain and the corresponding non-autonomous 
dynamical system is given by 



= A(t) 



-I 
I 



(24) 



where A(t) is a time-dependent strain amplitude. When A = const., the point (x, y) — (0, 0) is 
a hyperbolic saddle with a ID stable and unstable manifolds aligned with, respectively, e x and 
e y . When dA/dt j^— and A(t) > 0, it can be easily verified that -y(t) = is a trajectory of 



( |24[ ) in the extended phase space (x, y, t). Moreover, <y(t) is hyperbolic and has a 2D stable and 
unstable manifolds in the extended phase space which are spanned by, respectively, {e^et} 
and {e y ,et}. The fundamental solution matrix, X(t,to), of (241 is given by 



X(t,t ) 



p -4(Mo) 







(25) 



where A(t,t ) = f^A(r)dr. 

Note that the finite-time Lyapunov exponents, A1.2 (cf Definition A.l 1, for the flow associ- 
ated with ( 25 I are given by 



A T ' (x,y,t ) 



A(t + T, t ) 



2\T\ 



(26) 



and are independent on the spatial coordinates. Consequently, the FTLE field given by 
\T(x,y,to) = max[A 1 (a;, y, to), X 2 (x, y, to)] is spatially homogeneous and does not reveal any 
structure despite the fact that the stable and unstable manifolds of the hyperbolic trajectory 
■y(t) = are well defined. 



The Arnold cat map: The Arnold cat map, defined on the torus, is given by 



p n+1 = p n + q n (mod 1), (27) 
q n+1 = p n + 2q n (mod 1), (28) 

This dynamical system has a number of remarkable properties that are amenable to explicit 
analysis resulting from the linearity of the map and the doubly periodic boundary conditions. 
In particular, every trajectory can be shown to be hyperbolic and explicit expressions for its 
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Figure 5: Geometry of two material surfaces in the extended space (x, y, t) approximating the unstable manifold 
(red) and the stable manifold (blue) of the trivial solution, x(t) = 0, of the system For the chosen form of the 

amplitudes As, -4"gr ( c f P5|), the trivial solution is (infinite-time) hyperbolic on / = R but finite-time hyperbolic 
only on / = (— oo, —4.47] and / = (4.47, oo] (see text for definition of finite-time hyperbolicity on an interval). 



stable and unstable manifolds can be computed. The map can be shown to be ergodic, mixing, 
and to have the Bernoulli property, and each of these properties is present on the entire domain 
of the map. The proofs of these results are "well-known", but are often difficult to track down 
in the literature. [80| contains proofs, and also a guide to the original literature. However, for 
our purposes here we are only concerned with the Lyapunov exponents of trajectories of the 
cat map. These can be explicitly computed from the map and are found to be 

Ai, 2 = ±ln(3 + >/5)/2, (29) 

and they are the same for every trajectory. Therefore, we have a situation where, in some sense, 
the map is the "most chaotic possible" (i.e. it has the Bernoulli property) on its entire domain 
and every trajectory is hyperbolic (having one Lyapunov exponent with modulus greater than 
one and one Lyapunov exponent with modulus less than one) with stable and unstable manifolds 
that can be computed explicitly. Nevertheless, since the Lyapunov exponents of every trajectory 
are identical then contours of the FTLE fields are all identical, and this they reveal no LCS'aj 

Summary: We have shown two examples where the Lyapunov exponents can be explicitly 
computed for every trajectory. In each example the Lyapunov exponents were shown to be 
identical for every trajectory. Dynamically, these two examples could not be more different. 
The flow defined by a linear, time-dependent strain on the plain does not possess complex 
dynamics, even though (almost) every trajectory has a positive Lyapunov exponent. The 
Arnold cat map defined on the torus is extremely chaotic on its entire domain (and every 
trajectory also has a positive Lyapunov exponent). Clearly, complexity of trajectories is not 
sufficient for the FTLE field to reveal "structure". Rather, spatial heterogeneity is required, 
and this does not occur for linear flows, or flows exhibiting 'uniform' chaos, in the sense of 
identical Lyapunov exponents for (almost) every trajectory. 

5 This paper is concerned with an understanding of the role of manifolds and LCSs in fluid transport. Consequently 
we have been dealing with flows that are defined for continuous time. The Arnold cat map is a discrete time dynamical 
system. We have chosen it to illustrate a specific point because of its familiarity, and the ease for which its various 
properties can be explicitly computed. Nevertheless, the Arnold cat map dynamics can be realized in continuous 
time flows; see [TTJ [T2] EH for details. 



14 



Figure 6: (Left) Geometry (in the extended phase space (x, y,t)) of an unstable manifold of the trivial solution, 
x(t) = 0, in a flow generated by | |30[ )), (Right) Finite-time Lyapunov exponent fields, i.e. Xx(x,y,t) (cf |A.2| |, 
computed at three different times during the evolution t = 5 (top row), t = 7 (middle row), t = 13 (bottom row); for 
each of these times the FTLE fields were computed over two time intervals of different lengths T. The green lines 
denote the instantaneous geometry of the unstable manifold. When computed over sufficiently long time intervals, 
the ridges of the backward FTLE fields coincide with the unstable manifold. 



3.2.2 Strain-vortex-strain transition 

We consider here an example which is designed to illustrate the geometry and fate of finite- 
time stable and unstable manifolds of a finite-time hyperbolic trajectory during a flow transition 
associated with a loss and subsequent re-gain of finite-time hyperbolicity by this trajectory. We 
show here what kind of information about transport properties of such a flow can be obtained 
by analysing this transition using, respectively, the invariant manifold approach and the FTLE 
approach. 



Consider the following two-dimensional, non-autonomous dynamical system 



x = (A s {t)S(x)+A v (t)lT(x)^e- lM]2/s \ x 



G ST, t e M, 



(30) 



where 8 is a constant and the terms in the brackets represent a linear superposition (with 
time-dependent coefficients As{t) and Ajy(t)) of a straining field given by 



S(x) 



V 



and of a vector field with circular streamlines given by 



V(x) = 



(31) 



(32) 



Before proceeding to a discussion of concrete examples derived from (301, it is instructive 
to analyse the finite-time stability properties of the trivial solution, x(t) = 0. Some specific 
examples are discussed in the following subsection. 



Stability of the trivial solution, x(t) — 0. 



The linearisation of (30 1 about x(t) = is given by 



A(t)x = 


" -A s (t) 


-A v (t) ' 




X 










_ •%(*) 


A s (t) 




. y . 



(33) 
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\ T (x,y,t = 7) 



\ T (x,y,t - 5) 
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Fi gure 7'. Comparison of the ID sections (along (x,y = 0)) of the backward FTLE fields (gray-shaded) computed 
for the flow ( ] 30| t and discussed in figure [6] Three ID sections of the FTLE map computed at (a) t = 7 and (b) 
t = 5 over different integration intervals T. Note that the number of maxima and their location varies with T. 
In particular, the nature of the extremum at x = in the FTLE maps switches between minimum and maximum 
depending on T. The location of the finite-time unstable manifold of the (finite-time) hyperbolic trajectory x(t) = 
coincides, in this case, with the strongest maxima of the FTLE fields computed for T = 20. However, this fact can 
be only established once the finite-time unstable manifold is computed. 



Consider first a class of flows generated by (33 1 for which the coefficients, As(t), A~p(t) > 0, 
satisfy 



A~p[t) > Asit), for t 6 [t*,t**], -co < t* < t** < oo, 
■A-gr(t) < A s (t), for t G (oo,t*) U (t**,co). 



(34) 



In such a case, it can be shown that the trivial solution, x(t) = 0, t 6 M, has codimension-one 
unstable and stable manifold^] in the extended phase space {x,t). Consequently, it can be 
shown that the trivial solution is hyperbolic on 1R in the classical, infinite-time sense. How- 
ever, if we consider the finite-time stability properties of the trivial solution, some interesting 
issues arise. We note here that the theory of finite-time stability of non- autonomous dynami- 
cal systems is still an area of active research and, as a consequence, there exist, for example, 
at least two different ways of defining what is meant by finite-time hyperbolicity (cf Defini- 
tions A. 4 and A. 11 in the Appendix A). Although, it is currently not clear if these two notions 
are equivalent, or which one is more suitable for a given application, we show below that they 
predict essentially the same stability changes in the configuration considered here. 

As discussed briefly in the Appendix A, one approach to characterising finite-time stability 
properties of a given trajectory is via the notion of finite-time exponential dichotomy which 
is associated with a system linearised about this trajectory. While this notion of finite-time 
hyperbolicity seems more general and is very useful in more abstract considerations, it is often 
difficult to verify in practice. Nevertheless, provided that A~p(t) and As(t) are bounded and 

6 W e skip the proof here but the existence of the 'infinite-time' stable and unstable manifolds of the trivial solution 
of l |33| can be shown by using tech niqu es analogous to those used in 1541 cf §4]. The main difference here is the presence 
of the off-diagonal terms in A (cf |33k) which invalidates the contraction mapping argument when A~p(t) > As(t)- 
However, one can show the existence of a codimension-one manifold of trajectories on I = (—oo,t*) which converge 
to x = as t — > — oo, In the linear case of ( | 33 [ ) these solutions can be extended to / = (t* , oo) with the help of the 
fundamental solution matrix. Similar procedure can be used to show existence of trajectories of | |33| on / = (i**, oo) 
converging to x = as t — + oo, and then mapping them backwards using the fundamental solution matrix. 



1G 




-4 -2 2 4 -4 -2 2 4 -4 -2 2 4 



X X x 

Figure 8: A sequence of backward FTLE fields (grey/red shaded), Xx(x,y,ti), {ti\i<zi (cf |A.2| , for the flow ( |30| ) 
with As = 1 and A~^j(t) given by |46|. The FTLE fields are computed with |T| = 10. The flow undergoes a 
transition associated with the loss of finite-time hyperbolicity by the trivial solution. The dashed green lines denote 
the instantaneous geometry of a material curve which approximates the finite-time unstable manifold of x(t) = 
before the transition (After the transition the trivial solution does not have ft unstable manifold but this curve 
remains a material transport barrier in the flow.) Note that, when computed with a fixed T, the ridges of the FTLE 
field fade away during the evolution as the flow transitions into the 'non-hyperbolic' phase. 



sufficiently slowly varying on a time interval I, it can be shown [191 Propositions 1-2, p. 50] 
that the trivial solution is finite-time hyperbolic on any time interval J C I within which the 
real parts of the eigenvalues of the matrix A(t) in (33 1 are non-zero and have opposite signs. 



Conversely, it can also be shown [191 Proposition 2, p. 54] that a trajectory cannot be finite- 
time hyperbolic if the eigenvalues of A(t) are imaginary over a sufficiently long time interval 
(the slower the variation of the coefficient matrix, the longer time interval needed) . Since the 
eigenvalues of the matrix in (|33[) are given, at any t £ M, by 



(35) 



one can conclude that, if As and Ajy satisfy (34 1 and / = [t*,t**] is sufficiently long, the 
trivial solution is not finite-time hyperbolic on I. 

Another approach to characterising the stability properties originates from the so-called 
EPH-partition due to Haller (see the Appendix and [381 128| ). This criterion relies upon con- 
sidering the characteristics of the so-called rate of stra in tensor, S(t) (cf Definition A.7l, and 
the strain acceleration tensor, M(t) (cf Definition A. 9 1, derived for a flow linearised about the 
considered trajectory. In particular, a trajectory is said to be in a hyperbolic region of the 
phase space within a time interval / if the restriction of M{t) to the so called zero-strain set 
(cf Definition A. 8 1 is positive definite for all t 6 I. In the case of our example system ( 33 1, the 



rate of strain tensor 



S(t) = UA(t) + A(t) T ) 



-As(t) 
As(t) 



(36) 



is indefinite for any t 6 M, and the zero-strain set, defined as Z(t) = {a; £ R 2 : {x, S(t)x) = 0} 
is given by 



a e R 



(37) 
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Finally, the strain accelaration tensor is 



M(t) = S(t) + S(t)A(t) + A(t) T S(t) = 2 



A s {tf As(t)A v (t) 
A s (t)A v (t) A s (t) 2 



and its restriction to the zero-strain set yields 

{C,M(t)C) = a 2 As(t)(A 3 (t)-A^(t)), 

it, M(m + ) = a 2 As(t)(As(t) +A v (t)). 



(38) 



(39) 



(40) 



Consequently, the restriction of M(t) to Z(t) is positive definite provided that Asit)— A^(t) > 0. 
If the amplitudes As(t), A~^(t) satisfy (34 1, one can conclude, that the trivial solution leaves 
the hyperbolic region of the phase space at t* and is contained in the elliptic region (cf Def- 
inition A. 10 1 for t £ I = [t*,t**]. According to Definition A. 11 the trivial solution will not 
be finite-time hyperbolic on any time interval J £ IR such that J n I is sufficiently long (see 
Definition | A . 1 1 1 for more details). 

Note that both of these characteristics of finite-time hyperbolicity depend on the time 
interval considered and cannot be attributed to a point on a trajectory. Rather, whether or 
not a given trajectory is finite-time hyperbolic on a given interval, 7, depends on the relative 
length of subintervals of / within which the local dynamics has 'undesirable' properties. In 
what follows we will say that a trajectory 7 is not finite-time hyperbolic on an interval / if 
there exists interval(s) J such that J n I 7^ and 7 is not finite-time hyperbolic on J. Clearly, 
if a trajectory 7 is finite-time hyperbolic on / £ IR than it is finite-time hyperbolic on any 
J C I. 



Note also that if, instead of ( |34[ ), the amplitudes were chosen such that 

A*p(t) < As{t), for te(-oo,t*), 1 

A~p(t) > As(t), for t€(t*, 00), J 



(41) 



one can only identify 6 an unstable manifold, W u [x = 0], in the flow generated by (30 1. In such 
a case 5fte[<r(i)] = for any t £ [t*,oo] and the trivial solution is not finite-time hyperbolic on 
t £ [t* , 00) (i.e. x(t) = does not have the exponential dichotomy on [t* , 00)). Similarly, when 



Ay;(t) > As(t), for t £ (-00, t*) 
A~p(t) < As(t), for t e (t*,oa), 



(42) 



0] . The trivial solution is in this case finite-time 



one could only define a stable manifold W s [x 
hyperbolic on (t*, 00) but not on (— 00, t*). 

Note finally that if we restrict the system (30 1 to a bounded time interval, I = [t a , tb] C M, 
with t a > — 00, tb < 00, it is not possible to defin^] the stable and unstable manifolds of 
the trivial solution (in the classical, time-asymptotic sense) even if x(t) — is hyperbolic for 
the same system considered on I — JR. This situation is by far the most common one in 
applications, especially when dealing with experimentally measured or numerically generated 
flows. However, if x(t) = is finite-time hyperbolic on / (in the sense of Haller [15]), one can 
define (cf [28]) the following two flow-invariant sets: The t-Sbre of a finite-time stable set of 
x(t) = on / is given by 



WJ [x(t) = 0] (t) = jx £ R 2 : ||X(m, t)x\\ < 0, m £ /| , 
and the finite-time unstable set of x(t) = on / is defined, for t £ 7, as 



/" [x(t) =0]{t) = \xe IR 



dm 



\X(m,t)x\\ > 0, me I 



(43) 



(44) 



7 For systems defined on a finite time interval one can still consider non-unique extensions to I = IR by applying 
the Lyapunov- Perron approach to an extension of the flow from I = [a, b] to M as in |41II36I . Since such extensions 
can be accomplished in a non-unique way, the manifolds constructed in the extended system are unique up to an 
error C(e' cl - b ~ a ^),c> 0. 



18 



In contrast to the classical (time asymptotic) definition of stable and unstable manifolds, the 
finite-time counterparts, W" and Wf , have the dimension of the extended phase space (rather 
than a lower dimension) and their t- fibres are open sets in Ft 2 . In such a case, a common 
approach used in the invariant-manifold Lagrangian transport analysis is to choose (non-unique) 
segments of initial conditions of length a< 1, ilf a and ©£, containing the trivial solution of 
the linearised systenj^] and follow their forwards and backward time evolution. It can be shown 
(see Appendix[B| that, if properly chosen, the material segments are contained in, respectively, 
W" [x(t) — 0] and Wf [x(t) = 0] . Moreover, due to the the embedding property of finite-time 
stable and unstable manifolds (see |28l Theorem 37, p. 659]) the effect of the non-unique choice 
of the initial material segments diminishes with the length of the considered time interval /, 
provided that the considered trajectory is finite-time hyperbolic on / (see the Appendix [B] for 
more details). 



Examples of flows generated by the system ( 30 1 



In our comparison of the invariant manifold and the FTLE analysis of flows generated by 
(30 1, we first choose the amplitudes As and Ac in such a way that the flow is not finite-time 



hyperbolic on a bounded interval; this can be achieved, for example, by setting As(t) = 1 and 

Aff(t) = ^atan(10 - t 2 /2) + , (45) 

in which case t* « —4.47 and t** « 4.47 and the trivial solution is not finite-time hyperbolic 
on I — [t*,t**\. The results to computed for such a flow are discussed in figures 5)7 



In figure [6] we show the geometry of the numerically approximated unstable manifold of the 



trivial solution in the nonlinear flow ( 30 1 and compare these results with the backward FTLE 
fields (cf |A.2[ ) at three different times during the evolution t = 13 (top row), t — 7 (middle row) 
and t = 5 (bottom row). The unstable manifold was approximated by following an evolution 
of appropriately chosen initial material segment (cf Appendix |B|, using algorithms analogous 
to those described in [611 163] . Clearly, for sufficiently long integration intervals the ridges of 
the backward FTLE field coincide very well with the instantaneous geometry of the unstable 
manifold (dashed green), as can be seen in the panels computed with T = 20 at three different 
times (left column). Note, however, that for smaller values of T not only the ridges of the FTLE 
fields become less localised but their location changes as well. This effect is further highlighted 
in figure [7] where we show ID cross sections of the FTLE fields computed for different values 
of T. The non-uniqueness of the backward FTLE fields is a direct consequence of the fact that 
if one computes separation of nearby trajectories in non-autonomous flows, the outcome will 
depend, in general, on the starting time and the extent of the time interval over which such a 
diagnostic is evaluated. Therefore, in more complex flows it might be not always clear which 
length, T, of the integration interval is the most suitable one for describing the flow structure 
based on the FTLE fields. It is also worth noting that in complicated flows, possibly known 
on only for a finite time, the identification of Distinguished Hyperbolic Trajectories on a finite 
time interval and their stable and unstable manifolds is also not unique, although for different 



reasons (see [45] and the discussion following (431, ( 44 1 ) 



We finish this section with an example of a flow associated with ( 30 1 with As (t) = 1 and 



A-p(t) = - ( atan(10t) + tt/2 ), (46) 



which corresponds to the case (41 1 mentioned above with t* « 0. In figure [8] we consider a 
hypothetical situation of trying to record the time-dependent geometry of a transport barrier, 
given by the unstable manifold of x(t) — 0, using the backward FTLE fields. Note that, as 
discussed earlier, the trivial solution is not finite-time hyperbolic on any interval contained in 
I = (t*,oo) and, consequently, it does not have a finite-time unstable (or stable) manifolds 
on /. Assume that we choose a time interval length T which leads to well localised ridges in 
the backward FTLE fields during the initial period of evolution. In this case \T\ — 10 seems 
satisfactory for determination of the LCS before the transition. Nevertheless, it can be seen 
that the ridge localisation deteriorates in the FTLE fields, XT(x,y,ti), computed at an ordered 
sequence of 'observation' times {ti}i 6 Z with increasing tj. 



8 Note that, by construction, the trivial solution £(i) = of a system linearised about some trajectory r/(t) 
corresponds to this trajectory. 
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Figure 9: Backward FTLE field (cf |A,2[ | and invariant manifold structure for the double gyre flow (J47J) in the steady 
case, a) The most pronounced ridge of the FTLE map (see b) coincides with the hetcroclinic connection between 
the two hyperbolic fixed points pi and p2. Note that the FTLE field (b) possesses an inward spiraling ridge in each 
cell which does not correspond to an invariant manifold. Depending on the characteristics of the colour map, this 
spiralling structure can be suppressed. However, in a time-dependent case it is not immediately clear whether or not 
similar spiralling ridges correspond to transport barriers. 



3.2.3 Double gyre flow 

The double gyre flow is considered in the domain D — [0, 2] x [0, 1] and is given by 

u(x,y,t) — —nAsm(-Kf(x,t))cos(iry), 

df 



v(x, y, t) = ty A cos(nf(x, t)) sin(ny) 



dx ' 



(47) 
(48) 



where f(x, t) is chosen in such a way that /(0, t) — f(2, t) = 0. This flow is frequently used for 
illustrating the LCS (e.g. |75l 178] ) and it is instructive compare the LCS and stable/unstable 
manifolds of hyperbolic trajectories using this example. The (non-autonomous) dynamical 



system associated with ( 47 1 is simply given by 

x = u(x,y,t), 
V = v(x,y,t). 



(49) 



When the flow is steady, i.e. when df(x,t)/dt — 0, there are two hyperbolic stagnation 
points in the system (491 located at pi (x, y) — (1,0) andp2(:r,y) = (1, 1). The unstable mani- 
fold of the stagnation point pi coincides with the invariant boundary (a; = [0, 2], y = 0) and the 
stable manifold is located within the domain. Similarly, the stable manifold of p2 is contained 
in the flow-invariant boundary (x = [0, 2],y = 1) and its unstable manifold coincides with the 
heteroclinic connection between pi and p2. When f(x,t) = x the heteroclinic connection is 
given by (a; = 1, y = [0.1]). The steady situation is visualised in figure [9] where we overlap the 
forward FTLE field with the manifold structure associated with the hyperbolic fixed points pi,2- 
Unsurprisingly, the most pronounced FTLE ridge coincides with the heteroclinic connection 
discussed earlier. 

When df(x,t)/dt ^ 0, the paths of instantaneous stagnation points, pi(t), P2(t), are not 
system trajectories. They can, however, be used to compute two Distinguished Hyperbolic 
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a) 3L 1 fi(«,»,t = S9 b) \i 5 (x,y,t = 0) 




Figure 10: Backward (a,c) and forward (b,d) FTLE maps for the double gyre flow {47} at t = 0; computed over 
two different time intervals with lengths T = 15 (a,b) and T = 20 (c,d), At = 0.01. The parameters uj = 27r/10, 
e = 0.25 and A = 0.1 are chosen as in the online tutorial |75| . The dashed black curves denote the instantaneous 
geometry of the unstable manifold (a,c) and of the stable manifold (b,d) of the hyperbolic trajectories 72 (t) and 71 (t). 
For sufficiently long integration times, good agreement between the LCS (red) and the manifolds can be achieved. 
However, depending on T the FTLE map reveals ridges of different length and connectivity. Some most significant 
differences are marked by the black arrows. The correlation between the LCS and the invariant manifolds depends 
also on the integration method, the integration step At (see figure [ill. 



Trajectories (DHTs, cf Definition A. 6 Appendix A) which are contained in the flow-invariant, 
top and bottom boundaries. These DHTs can be computed using techniques described in 
[451 1631 161] . We stress again that the path of ISPs is just a convenient but not a necessary 
choice of the initial, frozen-time hyperbolic guess. In figure [10] we show examples of backward 
FTLE fields (a,c) and forward FTLE fields (b,d) computed at a fixed time, t = 0, over different 
lengths of the integration time interval T. We compare these results with the instantaneous 
geometry of the unstable manifold of a DHT 72 (t) (confined to the top invariant boundary) and 
of the stable manifold of a DHT 71 (t) (confined to the bottom boundary); the instantaneous 
snapshots of these manifolds are delineated by the dashed black curves. In the computations 
we used 

f(x,t) = a(t)x 2 +b(t)x, (50) 
a(t) = esinwt, (51) 
b(t) = 1 - 2esinu>t, (52) 

with u> = 27r/10, e — 0.25 and A — 0.1, which coincides with the choice used in [75]. Since 



the flow (471 with f(x,t) given by (50 1 is time-periodic, both the DHTs and their stable and 
unstable manifolds are well defined and unique. Moreover, since these manifolds are composed 
of the system trajectories, they represent material curves at any fixed time and are therefore 
barriers to transport. It can be seen that in this case the instantaneous geometry of the stable 
manifold of 71 (t) and the ridge of the forward FTLE map (i.e. the repelling LCS) are well 
correlated over long distances (in the arc length sense from the DHT). Similarly, the unstable 
manifold of 72 (t) and the attracting LCS associated with the backward FTLE map coincide 
provided that the FTLE field is computed over sufficiently long time interval. Recall that, 
as already mentioned earlier, at each 'observation' time the FTLE field, \T(x,y,t), depends 
on the integration parameter T. Thus, the arclength of the strongest ridges of the FTLE 
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Figure 11: Sensitivity of the FTLE field to the integration method. Forward FTLE computed at t = with T = 15 
for the flow (47]l using (a) 4th order Runge-Kutta and (b) forward Euler (used in the LCS MATLAB Kit [20] ); 
At = 0.1 in both computations. The fact that the results depends on the integration method and the time step used 
are hardly surprising. However, it is important to bear these effects in mind, especially when analysing experimental 
data (recorded on a discrete space-time grid) when one does not have the control over the choice of the discretisation 
(see also figure [l4] l. 



field and, more importantly, the location of these ridges varies with T. This can be seen in 



figure 10 c,d) which is computed for the same values of the flow parameters as in figure 10 a,b) 
but for T = ±20. Note, in particular, the changes in the FTLE fields occurring in the regions 
indicated by the arrows. 

Another interesting aspect related to the FTLE computations is the identification the LCS 
(i.e. the ridges of the FTLE fields) and their connectivity. The ridge extraction was described 
in [7H| and an example of the use of such a procedure can be found in [641 Figure 2]. However, 
it seems that such a ridge extraction is not commonly carried out. We note, for example, that 
the results discussed in [761 1771 175| and some results in |78j seem to be obtained not by ridge 
extraction but by appropriate 'thresholding' of the colour map used for shading the FTLE 
fields. In figure [12] we show a few examples of different shading of the same FTLE field which 
reveal a 'ridge landscape' of varying complexity with a number of disconnected ridges appearing 
(or disappearing), depending on the colour map threshold used. 

In summary, we observe a good correlation between the stable and unstable manifolds of 
DHTs and the ridge segments identified in the forward/backward FTLE fields in this flow. 
However, for a given FTLE field, the choice of the parameters T and the filtering applied 
to extract the LCS is rather subjective and can be ambiguous. This is of particular concern 
when analysing transport in time-dependent flows via the lobe dynamics. Lagrangian transport 
analysis via the lobe dynamics requires the ability to follow the evolution of lobes associated 
with tangles of stable and unstable manifolds of relevant hyperbolic trajectories. Any numerical 
technique for identifying these tangles will provide, at best, a good approximation of these 
structures. However, the minimum requirement for this kind of analysis is that the numerical 
method has to be capable of identifying evolution of the same, and sufficiently long, segments 
of such structures. If the structure of the relevant stable and unstable manifolds of DHTs is 
known, it is generally possible to adapt T and the colour map 'threshold' so that sufficiently 
long and connected LCSs are revealed. However, if the manifolds structure is not known a 
priori, this task may quickly become impossible. We also note that, while the methods based 
on computation of stable and unstable manifolds are capable of identifying and following in 
time long segments of hyperbolic structures, the necessary identification of the 'distinguished' 
hyperbolic trajectory is not always easy. Thus, it is likely that a synergetic approach, combining 
the use of FTLEs for identifying the possible locations of the DHTs with a subsequent manifold 
computation, may offer the right way forward. 

3.2.4 Time-dependent Hills' spherical vortex in the symmetry plane 

Consider now a class of velocity fields obtained by perturbing the well known steady solution 
of equations of an inviscid incompressible fluid flow given by the Hill's spherical vortex (see, 
for example, [7]). The Hill's vortex flow, H, is then perturbed by a time-dependent strain, S, 
so that the corresponding dynamical system is given by 

x = H(x,y,z) + S(x,y,z,t). (53) 
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Figure 12: When shading the FTLE fields one has to make a choice of a filtering threshold in order to reveal the 
ridges approximating the LCS (backward FTLE field at t = computed for T = — 20 shown). Different choices of the 
colour mapping, which serves here a height filter, may reveal or suppress disconnected segments of LCS. This effect 
combined with the non-uniqueness of the FTLE maps (i.e. the one-parameter family {Xt(x, y, t)} t +Tel) makes it 
difficult to identify long segments of the LCS which are necessary for transport analysis via the lobe dynamics. 



The components of the steady Hill's vortex in Cartesian coordinates are 

H x ~ (itr sin 6 + uq cos O) cos <E>, 
H y = (u r sin B + ue cos O) sin $, 
H z = (iir cos O — sin 0) , 



(54) 



where r = ^ x 2 + y 2 + z 2 , B = acos(z/r), $ = a,cos(x / ^/x 2 ^{-y 2 ) and, assuming that a denotes 
the radius of the vortex, the velocity components in the spherical coordinates are 



U(l - a 3 /r 3 )cos6 if r^a, 
-|(7(l-r 2 /a 2 )cos6 if r < a, 

-[/(l + a 3 /(2r 3 ))sin6 if r > a, 

(J(7(l - 2r 2 /a 2 )sin6; if r < a. 
This unperturbed (steady) Hill's vortex flow has two hyperbolic stagnation points 

fti = (0,0, -a) T , h 2 = (0,0, a) T , 



«e 



(55) 



(56) 



(57) 



which are located on the (flow-invariant) axis of symmetry e z of the vortex. The fixed point 
hi has a two-dimensional unstable manifold in lR Li (ID in any symmetry plane containing e z ), 
and the fixed point I12 has a two-dimensional stable manifold in M 3 . 
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Figure 13: Schematic representation of a three-dimensional flow used in computations of invariant manifolds 
and FTLE fields in 1 3.2.4 The steady Hill's spherical vortex (a), sketched in a symmetry plane, is perturbed 
by a time-dependent strain (b). One of the principal axes of the straining flow is aligned with the axis of 
symmetry of the Hill's vortex, e z . The amplitude of the strain changes with time as shown in c). 



The perturbing, time-dependent straining 


flow 


is given by 




a(t) 










X 


S = A(t) ■ 





(3(t) 







y 










7(*) . 




z 



(58) 



where A(t) is a time-dependent amplitude, the strain rates are normalised so that max(a, j3, 7) = 1 
and they satisfy a + + 7 = 0. 

When < A <C 1 the fixed points hi and hi no longer exists but they are perturbed to 
two hyperbolic trajectories, 71 (t) and 72 (t), which possess, respectively, a 3D unstable and 3D 
stable manifolds in the extended phase space spanned by {e x ,e y ,e z ,et} ■ In other words, at 
any fixed time instant the unstable manifold of 71 (t) and the stable manifold of 72 (t) are given 
by surfaces embedded in M . 



As long as one of the axes the perturbing straining flow ( 58 1 is aligned with the symmetry 



axis of the Hill's vortex, the flow (53 1 remains axisymmetric. Consequently, every plane con- 



taining e z is invariant with respect to the flow ( 53 1, with H given by ( 54 1 and S given by ( 58 1 



We therefore restrict the analysis to one such symmetry plane, namely (x = 0, y,z), in which 
the instantaneous geometry of the considered invariant manifolds is given by curves. 

The hyperbolic trajectories, 71 (t) and 72 (£), which are confined to the symmetry axis, 
e z , can be computed using the same algorithms (cf [451 50 ) as used in the previous examples. 
Their stable and unstable manifolds are computed as in the previous examples using techniques 
described in |63II61| with the initial 'seed' for these computations chooses in the way described 
in Appendix B. In figure [14] we compare the instantaneous geometry of the unstable manifold 
of 71 with the corresponding backward FTLE field, both computed in the symmetry plane for 



the flow associated with ( 53 \ with the perturbing strain amplitude given by 



A(t) = (0.05 + 0.3sin(2.33i))e 



-(t-l) 2 /(3.5) 2 



(59) 



The strain rates are chosen as a — (5 — —0.5, 7=1. The conclusions one may draw from these 
computations are similar as those drawn from the previous examples. Provided that the FTLE 
fields are computed with sufficient care the overall agreement between the ridges of the FTLE 
field (red) and the unstable manifold is rather striking (cf figure 



14 



a)). As in the previous 
examples, the critical parameters for an accurate manifold computation are the maximum and 
minimum curvature cut-off parameters and accurate integration routine. 

However, we intend to use this flow geometry to alert the reader to the potential problems 
which are particularly likely to appear when analysing experimentally measured flow fields or 
velocities obtained from numerical PDE solvers. 

In order for the FTLE computations to be reliable, one needs to make sure that the com- 
putational grid is sufficiently refined to reveal the desired details and, most importantly, that 
the integration routine is chosen appropriate for the chosen integration time step. Obviously, 
in the case of analytically defined flow fields, as the ones we are dealing with here, the choice 



24 




Fi gure 14: Comparison of the backward FTLE fields, Xx{x^y y t) (cf |A.2[ >, computed with T = —15 in a symmetry 
plane for the axisymmetric, time-dependent, perturbed Hill's vortex flow (|53[l a t t = 0, and the instantaneous 
geometry of the unstable manifold of the hyperbolic trajectory 7i(t) (see |3.2,4| l. The top row shows the FTLE 
fields computed using 4th order Runge-Kutta with a) At = 0.01 and b) At = 0.1. The bottom row shows analogous 
computations performed using the forward Euler method (as in |20p with c) At = 0.01 and d) At = 0.1. Provided 
that an appropriate method is used for the integration of trajectories (i.e. not the forward Euler) a good agreement 
can be achieved (cf (a)) between the LCS computations and the invariant manifolds computations. 
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of the integration time step is not a serious constraint. However, in the case of discrete data 
sets (numerical or experimental) the time-discretisation of the data set imposes limitations on 
At, requiring a trade-off between the time step chosen and the temporal data interpolation. 
In order to highlight, the kind of problems one might encounter in such a situation we show, 



in figures 14 c,d) results of the FTLE computations for the same flow as in figures 14 a,b) 
but using the first-order accurate forward Euler integration method. This method is in fact 
implemented in the LCS MATLAB Kit mentioned earlier [20] which is combined with linear 
spatial interpolation of the discrete flow data required by the code. Note, in particular the 
erroneous structures in figure [l4^ d) which emerge in the FTLE fields computed using the for- 
ward Euler integration method with At = 0.1. The main danger here is associated with the 
main advantage of the FTLE computations. Namely, it is straightforward to develop a basic 
algorithm computing FTLE fields which will generate reasonably looking output. 



3.2.5 Boundary layer separation on a non-slip boundary 

The technique of invariant manifolds and lobe dynamics for finite-time, aperiodically time- 
dependent velocity fields has not been extensively developed. An important area of application 
in this setting is separation from a non-slip boundary. In this setting Haller and co-workers 
have developed a comprehensive theory based on the FTLE and LCS approach [841 1401 
1811 182] - Related earlier work using non-hyperbolic separation points and manifolds can be 
found in [791 1271 1871 133] . Nevertheless, there has been extensive work in the mathematics 
literature on non-hyperbolic trajectories and their stable and unstable manifolds, e.g. [651 1161 
1301 1171 1151 [5j HOI 1441 [6] . This work should serve as an excellent foundation for developing a 
theory of 'distinguished saddle-points' and their stable and unstable manifolds in finite time, 
aperiodically time-dependent velocity fields. Finally, we note that the algorithm for computing 
time-dependent invariant manifolds described in 63, 61 does not require a hyperbolic trajectory 
as a starting point. Rather, it requires an appropriate 'seed' from which the material curve, 
approximating an invariant manifold is 'grown' according to the numerically integrated vector 
field. Depending on the choice of the 'seed', the obtained results may, or may not, be relevant 
for transport considerations. Instead of selecting the location of some distinguished hyperbolic 
trajectory as the 'seed', one could choose the instantaneous location of a non-hyperbolic saddle 
point. However, this situation has yet to be developed. 



3.2.6 Eddy-pair system 

In this example we focus on a flow exhibiting a transition between a configuration characterised 
by a single Lagrangian eddy and an eddy pair. This is a simplified version of the kinematic 
model of the front-eddy system introduced earlier in |13| . As in all other examples in this 
section, our main objective is to establish how well the LCS, represented by ridges of the 
FTLE fields, correlate with invariant stable and unstable manifolds of DHTs in aperiodically 
time-dependent flows. 

The flow considered here is chosen in such a way that it undergoes a transition from a 
single Lagrangian eddy configuration to an eddy-pair configuration. The streamfunction of the 
kinematic model we use in our analysis is given by 

il> = M(t)+ ^ 1 (t) e -((--iW) 2 +(!'-^W) 2 )/*?( t ) + ^ 2 (f) e -((— i( t » 2 +^-fi( t )) 2 )/ 5 2 2 W ) ( 60 ) 

where 

M (t) = £(*)- a(t) (x cos Q(t) - y sin 9{t)) 2 + (3{t) (as sin 6(t) + y cos 9{t)) , (61) 

and C = — 1, a = 0.08, {3—1,0 = — it/A. The second and the third term in ( |60[ l give rise, for 
appropriate values of the amplitudes Ai and A2, to the appearance of closed contours in the 
instantaneous streamline patterns. We refer to such patterns Eulerian eddies. We choose here 
Ai = 10, Si = 4, xi = 2/1 = 4, X2 — y2 — 0, S2 = 0.9 and the time-dependent amplitude of the 
second Eulerian eddy as 

A 2 (t) = -2/7r(atan(t - 1) - atan(-9)). (62) 

With the above choice of the amplitudes Ai and A2 the flow is aperiodically time dependent 
and asymptotically steady, so that the two DHTs approach the location of the single fixed point 
in the system for t — > —00 and two fixed points at t — > 00. The DHTs are again computed using 
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the MATLAB implementation of the techniques introduced in |45l 149) . and their manifolds are 
computed using the ideas based on [631 161| . 



Figure 15 shows the backward and forward FTLE fields (yellow/red shades, see Defini- 
tion A. 2 1 computed for the flow (601 at t — 0, and stable (blue and cyan) unstable (magenta) 
manifolds of the two DHTs present in the flow. The location of this 'observation' time relative 
to the geometry of the DHTs is shown in the top left panel. The top right panel shows the 
backward FTLE map computed with T = 25 and an unstable manifold (of the two DHTs, they 
are extremely close if not identical). Clearly, the attracting LCS, corresponding to the ridge 
of the backward FTLE field, and the unstable manifolds of the two DHTs correlate very well 
over long arclength distances from the DHTs (black dots). The bottom panel shows a compar- 
ison between the stable manifolds (blue and cyan) of the DHTs and the forward FTLE map, 
showing a good agreement. Note also the spiral structure in the forward FTLE map (bottom) 
which is visible inside the small eddy. When computed over long time intervals the length and 



definition of the extracted ridges might increase (see, however, S 3 . 1 1 but the method starts 
detecting 'premonitions'/'ghosts' of the future/past phase space geometry. Note also that the 
significant inward curl of the LCSs inside the large eddy in both forward and backward FTLE 
fields which does not correspond to the manifold geometry. 

We show two more snapshots of the instantaneous geometry of the FTLE fields and the 



stable and unstable manifolds of the DHTs at t = —4 (figure 161, and at t = —8 (figure 17 1. 
In all cases the agreement between the dominant FTLE ridges and the corresponding stable 
or unstable manifolds of the DHTs is good, provided that the FTLE fields are computed for 
sufficiently large T. 

3.2.7 Eddy-quadrupole system 

In this final example we focus on an incompressible flow characterised by the following stream- 
function 

i/)(x,y,t)= (xy(o(t)-x 2 ) - axy 3 + f3xy 5 ^je-^ +v ^' 5 \ (63) 

where a(t) is some function of time and a, /3, 8 are constants. The dynamical system associated 
with the flow is given by 

x = &ip/dy, y = -&tj)/dx, (x,y)&R 2 , t <E R. (64) 

We will choose here a particular form of time-dependence which will induce a symmetric tran- 
sition of the flow associated with (|64[) from a four-eddy configuration to an eight-eddy config- 



uration (see figure 181. We use this setting to illustrate two issues affecting, respectively, the 
invariant manifold computations and the FTLE computations. Due to the type of transition 
considered here we are not able to identify DHTs throughout the time interval considered. The 
problem affecting the FTLE computations stems again from their non-uniqueness and the fact 
that, in this case, the FTLE fields computed for longer integration times show less pronounced 
ridges, detecting ghosts of pre-transition flow characteristics. 



After a bit of algebra, one may notice that ( 64 1 has two invariant lines given by y = and 



x — 0. Alternatively, in the extended phase space one may identify two invariant planes 

I x = {(x,y,t) £ R 2 x R: y = 0}, l y = {(x,y,t) <E R 2 x R : x = 0}. (65) 
Note further that the dynamics in I x is given by 

y = Q, x = x(a(t)-x 2 )e~ xi/si = x(a(t) - x 2 ) + 0(x 5 ), (66) 
and the dynamics in X y is given by 

i = 0, y = -y{a(t)- ay 2 )e- vi/si = -y(a(t)-ay 2 )+0(y 5 ). (67) 

Clearly, we already analysed this type of one-dimensional dynamics in §3.1| In this example 
we will only consider the time-dependence that corresponds to Scenario II discussed in §3.1| 
i.e. we choose 

a(t) = 2(atan(10t) + vr/2 - l) ; (68) 
the remaining parameters in (|63[) are a = 1/3, (3 = 0.008/5 and 8 — 5. 



27 




b) 



Figure 15: FTLE fields (yellow and red shades, c fjA.2| and stable/unstable manifolds of two DHTs (black) computed 
in the flo w )|60| at t = (see (a)) during a transition between the singe-eddy and eddy-pair configuration (see the 
figures [To] 1 1 7| for the geometry at earlier times), (b) backward FTLE field, computed with T = — 25, superimposed 
with the unstable manifolds (dashed black) of the two DHTs (black dots); the LCS are delineated by the red ridges of 
the FTLE map and were enhanced by appropriate filtering of the colour map. (c) the forward FTLE field (yellow/red 
shades), computed with T = 25, superimposed with the stable manifolds (cyan/blue) of the two DHTs (black dots). 
The manifold segments inside the black rectangle were removed in order to reveal the LCS underneath. When 
computed over sufficiently long time intervals, the length and definition of the strongest ridges (LCS, red) of the 
FTLE maps generally increases (see, however, |3.1| l but the method starts detecting 'premonitions'/'ghosts' of the 
future/past phase space geometry. Note, in particular, the spiral structure inside the small eddy visible in the forward 
FTLE map (c). Note also a significant inward curl of the weaker ridges of the forward and backward FTLE fields 
inside the large eddy which are not associated with the instantaneous geometry of the invariant stable/unstable 
manifolds. 
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Fi gure 16; Comparison between the forward FTLE map (yellow/red shades) and stable (blue, cyan) manifolds of 
two DHTs (black dots) computed in the flow | |60[ l at t = —4. The FTLE field was computed with T = 25 and the 
LCS, represented by the red ridges, were 'extracted' by filtering the colour map. The manifolds inside the black 
rectangle were removed in order to reveal the LCS underneath. See figures [TB"] and [l7| for the geometry at other times 
during the transition. 
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Fi gurG 17: The forward FTLE map (yellow/red shades) superimposed with the stable (blue, cyan) manifolds of two 
DHTs (black dots) computed in the flow ([60]) at t = -8. The FTLE field was computed with T = 25 and the LCS, 
represented by the red ridges, were 'extracted' by filtering the colour map. The manifolds inside the black rectangle 
were removed in order to reveal the LCS underneath. See figures [15] and [TB] for the geometry at other times during 
the transition. 
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Figure 18: (centre) Dynamics in the two invariant planes, X x = {(x,y,t) e R 2 X R : y = 0} and 
X y = {(x,y,t) £ R 2 X R : x = 0}, in the extended phase space of the flow associated with (I64J) with a(t) 
given by j68l. The flow undergoes a transition associated with changes in finite-time stability properties of the 
trivial solution x(t) = (see text). The Distinguished Hyperbolic Trajectories are marked by thick black lines and 
paths of instantaneous stagnation points (ISPs) are marked by dashed green lines (and by green dots in (a-c)). The 
dynamics in the invariant plane I x corresponds to Scenario II in ^3.1 1 (a-c) Instantaneous streamline patterns in the 
flow associated with (|64b at three different times. 




Figure 19: Backward FTLE fields, X T (x, y, t) cf |A.2[ computed for the system |64| l with a(t) given by ((38| at t = 5 
and different integration time lengths a) |T| = 1, b) |T| = 3, c) |T| = 5, d) | T| = 10. ID cross sections of these fields 
along (x,y = 2) are shown in the central panel. The flow associated with (164b undergoes a transition which results 
in an emergence of four new eddies which are present in both the Eulerian and Lagrangian frameworks. Contrary 
to common intuition, the location of the strongest ridges in the FTLE fields varies with T and the overall strength 
of the ridges diminishes with T. This phenomenon is a direct consequence of the transition. Note also that the 
strongest ridge in (d), located at x = 0, is a 'ghost' of the the dominant repelling structure before the transition. 
The stable manifolds (blue) and unstable manifolds (red) of four hyperbolic trajectories involved in this process are 
shown in (e). 
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With a(t) given by (681 so that a(t*) = at t* fa 0.0642, one can easily see that within 



the plane I x the trivial solution, x(t) = 0, of (64 1 is pullback attracting (cf (141 and [55]) on 



/ = (— oo, t*) and that it is repelling (in the sense of (m)) on I — (t*,oo). If we consider the 
dynamics within the invariant plane I y , the trivial solution is repelling on I = (— oo,£*) and 



it is forwards attracting (cf ( 12 l and [55]) on I = (£*, oo). Consequently, while x(t) = is not 
hyperbolic on M (in the traditional, inifinite-time sense) it is certainly finite-time hyperbolic 
on any time interval which does not contain t* . Moreover, while for J C (— oo,t*) any finite- 



time unstable manifold of x(t) = 0, i.e. Wj[x(t) = 0] (cf Appendix |B| and { 3.2.2 l, contains 
a subset of T y , for any J C (i*,oo) the unstable manifold Wj[x(t) = 0] contains a subset of 
I x . The converse is true for the finite-time stable manifolds, Wj[x(i) = 0], for, respectively, 
J C (-oo,t*) and J C (t*,oo). 

Similarly to the one-dimensional dynamics considered in §3.1| the changes in stability prop- 
erties of the trivial solution are accompanied by a transition in the Lagrangian flow structure, 
which is associated with changes in the geometry of certain distinguished, hyperbolic trajecto- 
ries. Due to the presence of higher order terms in (|66|) and (|67|) we cannot compute the distin- 



guished trajectories in a way analogous to (10 1. However, one can resort here to the iterative 



algorithm (cf Appendix A and [451 149| ) as in most other cases discussed in this work. Recall 
that, as it was shown in [45, 45], if the iterative algorithm converges, it returns a hyperbolic 



trajectory. Such a trajectory is branded 'distinguished' if it is also bounded (cf Definition A. 6 1 
on the considered time-intervaQ Since we are concerned in this example with a system which 
is defined on / = IR and asymptotically autonomous due to the form of ( |68[ ), the obvious 
candidates for the location of the DHTs for t — ► ±oo are given by the hyperbolic stagnation 



points of the autonomous dynamical systems given respectively by (64 1 with a* = lim a(t). 

t — >±oo 



One would then expect the existence of five DHTs after the transition and one DHT before the 
transition. Since the DHT are, of course, trajectories, they cannot bifurcate. Consequently, 
all of the five trajectories which would be branded 'DHTs' after the transition must exist in 
the flow before the transition. For a given time interval I C M the finite-time DHTs can be 
located using the iterative algorithm provi ded that one can choose an initial guess, given by 



a frozen-time hyperbolic (cf Definition A. 5 1 path, which is C 1 and lies sufficiently close to the 
sought DHT (cf (821 in the Appendix A). Often, a useful initial guess can be constructed from 
the paths of instantaneous stagnation points which are frozen-time hyperbolic. This strategy 
is also useful here for finding two DHTs contained in the invariant plane I y . However, due 
to the nature of the dynamics in I x (which is identical with that considered in Scenario II of 



6 3.1 1 we are unable to construct a guess on intervals containing t* which would lie sufficiently 
close to the DHT. Identification of DHTs on intervals contained in (t* , oo) does not pose such 
difficulties but the outcome depends on the time interval chosen, i.e. the iterative algorithm 
converges onto different hyperbolic trajectories depending on the considered time interval (see 
the footnote). We compare the stable manifolds of the identified DHTs with ridges of FTLE 
fields computed for this flow in figure [T9[ 



When attempting to characterise the flow associated with ( 64 1 and ( 68 \ using the FTLE 
fields, one can, as in the previous examples, identify the one-parameter family of FTLE fields, 
{\t(x, y, *)} TgJR ! which are computed over different integration time intervals. Despite this 
non-uniqueness of the FTLE diagnostic, in most examples presented so far one could obtain 
good agreement between the invariant manifold calculations and the LCS obtained from At for 
sufficiently large T. In this case, however, the situation is rather different and in many ways 
analogous to the one-dimenisonal configuration discussed in |3.1| In figure [19] we show results 



of backward FTLE computations at t = 5 for the flow associated with ( 64 1 and o(t). The panels 
(a-d) show results of computations over four different lengths of the integration time interval 
(a) \T\ = 1, (b) |T| = 3, (c) \T\ = 5, (d) \T\ = 10; the central panel show ID cross-sections of 
these fields at (x, y — 2). Three issues affecting the ridges of the shown FTLE fields are worth 
noting: (i) the geometry of the ridges (i.e. the LCS) and their connectivity changes with T, 
(ii) the relative and absolute strength of the ridges diminishes with T, (iii) for sufficiently long 
(backward) integration times the strongest ridge in the FTLE field corresponds to a 'ghost' of 
the pre-transition flow structure (see y — in (d)). Consequently, it is rather difficult to obtain 
a coherent picture of the flow structure based on the family of FTLE's, {\t(x, y, *)} T6JR > m 
this case. 

9 Note that on a finite time interval this notion is non-unique since any trajectory of a smooth vector field is 
bounded on a bounded time interval. However, the ambiguities due to the non-uniqueness are, in general, only 
non-negligible near the end points of the time interval; cf |49l I45| . 
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4 Conclusions 



In this paper we have considered issues associated with the characterisation of the infinite time 
notion of hyperbolicity for aperiodically time dependent vector fields that are only known on 
a finite time interval by exploring concepts of finite-time hyperbolic trajectories, their finite 
time stable and unstable manifolds, as well as (one-parameter) families of finite-time Lyapunov 
exponent (FTLE) fields and associated Lagrangian coherent structures. Our approach has been 
to consider a collection of diverse examples where explicit phenomena can be exhibited and 
controlled. 

In Section |3. 1 1 we considered a one-dimensional vector field where the aperiodic time- 
dependence was specified in three distinct ways. This enabled us to probe the phenomenon of 
flow transitions and show how they may give rise to ambiguities in the effort to determine flow 
barriers from non-unique FTLE fields. Similarly, we used this configuration to illustrate issues 
associated with the lack of a unique, locally distinguished hyperbolic trajectories organising 
the structure of the flow. In Section |3.2.f | we considered two essentially dynamically opposite 
examples where the Lyapunov exponents of every trajectory could be determined analytically. 
In each example all Lyapunov exponents were identical, hence the FTLE fields did not give rise 
to LCSs. This highlighted the point that the emergence of LCSs is a consequence of spatial 
heterogeneity in the FTLE field and not just due to a rapid separation of nearby trajectories. In 
Section [3. 2. 2| we considered a velocity field exhibiting the strain- vortex-strain transition. This 
example illustrated some crucial issues associated with attempts to understand the nature of 
transport barriers in a transitioning flow in the finite-time setting. In this case, depending on the 
length and location of the 'observation window', different diagnostics could be obtained when 
employing the invariant manifold approach as compared with FTLE fields. In Section |3.2.3| 
we considered a double-gyre flow which has become a common benchmark flow in the LCS 
literature. We used this example to show that essentially the same information about the flow 
structure can be obtained from both techniques provided sufficient care is taken. This seems to 
be a common situation in flows which do not undergo transitions. We also illustrated there the 
sensitivity of the results to the order of the integrator used in computation of the trajectories 
as well as the importance of the cut-off level for the filtering procedure used in extracting LCSs. 
These conclusions and, in particular, the need for accurate trajectory integration, is further 
stressed in sect ion |5 . 2 . 4 1 where we considered an axisymmetric, time-dependent perturbation of 
the Hills spherical vortex. This flow serves as a good illustration of how inaccurate integration 
of flow trajectories can lead to plausible yet incorrect FTLE fields. The two closing examples, 
considered in sections |3.2.6| and |3.2.7[ were linked to the one-dimensional examples discussed 
in §3.1| The kinematic model of an 'eddy-pair system', discussed in §3.2.6[ is a common feature 
in geophysical flows and both the invariant manifold and the FTLE methods yield correlated 
diagnostics of the flow structure in this case. The 'eddy-quadrupole' system, discussed in §3.2.7| 
further highlights the problems that might arise when trying to select the most suitable FTLE 
field from the family parameterised by the integration time length. In particular this example 
illustrates the ambiguities one may encounter when attempting to increase the length of the 
integration time interval in order to obtain longer and more pronounced ridges in the FTLE 
field. Finally, we collect a number of technical details on finite-time hyperbolicity and its use in 
understanding fluid transport as well as a detailed discussion of an important technical detail 
concerning the choice of the initial material segment for the computation of finite time stable 
and unstable manifolds of finite time hyperbolic trajectories. 

The phenomena discovered and analysed in our examples point the way to a variety of 
directions for rigorous mathematical research in this rapidly developing, and important, new 
area of dynamical systems theory. 

Acknowledgements 

The authors acknowledge financial support from ONR Grant No. N000I4-0I-I-0769, and the 
stimulating environment of the NSF sponsored Institute for Mathematics and its Applications 
(IMA) at the University of Minnesota, where the work on this manuscript was begun. 



32 



A Some important definitions 



In order to make the discussion presented in this paper relatively self-contained, we recapitulate 
here some fundamental notions and definitions which are important for the analysis presented 
in the preceding sections. All of the material included in this section can be found in existing 
literature and we provide references, which are not exhaustive, to some relevant material. 

Consider a velocity field v : M n x I — > M n defined over a time interval I — [ti,tf] C JR and 
a system of ODE's 

x = v(x,i), iefi", tel. (69) 
The curves, f(t) : I M n that satisfy p9| , i.e. -y(t) = v("f(t),t), are referred to as M n - 
embedded trajectories of the non-autonomous dynamical system associated with v (other em- 
beddings are possible; for example, 7(f) : / — > ]R n x I, but we do not require such notions 
here) . 



which is not necessarily a trajectory of (691, in the following form 



Consider now a linearisation of the system (69 1 about a C r (r ^ 1) path, x(t) : I — > R n 



e = act>(3(t),t)e, (70) 

where £(t) = x(t) — x(t) and 9 x w(x(t), t) is the Jacobian of v(x, t) evaluated at x — x(t). We let 
X(t, ti) denote the fundamental solution matrix of (|70[), i.e. it is the map X(t, £;)(•) : Mt n — > iR n 



which is linear in t, and Lipschitz in t. Moreover, if £(t,£i,ti) is a solution of (70 1, then 
£(t,$i,ti) = X(t,ti)£i, and X(t,a)X(s,U) = X(t,tj). 

Consider first a situation when x is a trajectory of (|69[), i.e x = 7(f) and -y(t) = v{-i{t),i). 



Then, ( 70 1 describes the dynamics in the neighbourhood of the trajectory y(t) in the frame 
moving at speed 7. Thus, if <5t i denotes the perturbation of f(t) at t = ti, we find that it 
evolves according to 

||f(t)|| = y/(X.(t,t i )8t t ,X.(t,U)6t i ) = ^(8 u ,TL(t,U) T X.{t,U)8 u ), (71) 

where A = X(t, ti) T X(t, t,:) is commonly referred to as the finite-time Cauchy-Green tensor. 
Since A is real and symmetric, it can be diagonalised in an orthogonal basis of eigenvectors 
which denote the principal axes of growth of the infinitesimal perturbation. It then follows 
that the tensor 

M = (X(t,t l ) T X(f,t ! )) 1/2(t ~' l) , (72) 
is also diagonalisable in the same orthogonal basis. 

Definition A.l (Finite-time Lyapunov exponents, X T (x, t)). The logarithms of the eigenvalues 
of M and are called the finite-time Lyapunov exponents computed at time t over the time 
interval T. If T > 0, Xj<(x,t) is called the i-th forward finite-time Lyapunov exponent. If 
T < 0, \ l T (x,t) is called the i-th backward finite-time Lyapunov exponent. 

For more details regarding properties of Lyapunov exponents the reader is referred to ( [511 
1561 157] ). and for description of algorithms allowing their computation see, for example, ( |23l 

Definition A. 2 (Finite-time Lyapunov exponent field, Xt(x, t)). Assume that 

\ 1 T (x,t),\ 2 T (x,t),...,\^(x,t), (73) 



represent the finite-time Lyapunov exponents computed for a trajectory of (69 1 passing through 
x G M n at t. The scalar field 

X T (x,i) = max[\T(x,i), \r(x,t), . . . , \T(x,t)], (74) 

is called the finite-time Lyapunov exponent field at time t computed over a time interval of 
length T. If T > 0, it is called a forward FTLE field and if T < 0, it is called a backward FTLE 
field. 



Definition A. 3 (Finite-time Exponential Dichotomy). We say that the linear equation (701 
has an exponential dichotomy on the finite time interval / if there exists a (constant) projection 
operator P € ]R nxn , P 2 = P, and positive constants K, L, a, (3 such that: 

|X(t,t i )PX" 1 (s,t i )| < Ke~ a(t ~ s \ for t> s, t,s€l, 
|X(f,^)(Id-P)X _1 (s,ti)| < Le~ l3{s - t \ for s > t, t,s€l. (75) 
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For more details see, for example, 19, 43 . The notion of a generalised exponential dichotomy, 
where P does not have to be constant, is discussed for example, in [85]. Numerical methods 
for calculating the constants K, L, a, and j3 are given in [23] . 

Using the notion of exponential dichotomy, we can provide one possible definition of finite- 
time hyperbolicity. 

Definition A. 4 (Finite-time Hyperbolicity). We say that the path x(t) : I — > M n is finite-time 
hyperbolic on the interval I if the equation (701 has exponential dichotomy on I. Furthermore, 
if j(t) is a trajectory of the system (691 and we let x = 7(t), then 7 is called a finite-time 



hyperbolic trajectory on the interval / if the equation 

e = a.t;(7(t),t)6 



(76) 



has exponential dichotomy on /. 

Remark: In the limit tj — > — 00, tf — > 00 and x = 7(t), the above definition becomes equivalent 
to the standard notion of a hyperbolic trajectory. 

Roughly speaking, finite-time hyperbolicity implies that there exists a fc-dimensional (k ^ n) 
subspace in IR n of solutions tending to zero exponentially in forward time, and a (n — k)- 
dimensional subspace of solutions approaching 7(4) at an exponential rate in backward time; 
no assumptions are made the about the fate of these neighbouring trajectories beyond I even 
if the velocity field v(x,t) is known outside this interval. 

Definition A. 5 (Frozen-time Hyperbolicity). We say that the path x(t) : I —* FT 1 is frozen- 
time hyperbolic on the finite interval I if the eigenvalues of the Jacobian, d x v(x (t), t), in the 
linearised equation ( 70 1 have non-zero real parts for any fixed t 6 I. 

Remark: It can be shown, using results discussed in [49], that a path which is frozen-time 
hyperbolic is also finite-time hyperbolic (but not vice versa). 



Using definitions |A.3| and | A.4[ we can finally identify a special type of a hyperbolic trajec- 
tory which, if present in a considered flow, plays an important role in Lagrangian transport 
considerations. 



Let x (t) be a finite-time hyperbolic path and consider the nonlinear equation ( 69 1 in a frame 
'moving' with x by setting x(t) — y(t) + x(t). The transformed equation can be written as 



where 



y = A(t)y + f(y,t), y£]R n ,teI. 

A(i) = d*v(x(t),i), 
f (V, t) = v(y(t) + x(t)) - d x v(x(t) y t)y(t) - x(t). 



(77) 

(78) 
(79) 



Since we assumed that x(t) is finite-time hyperbolic (see Definition A. 2), we can associate 
the particular solution of ( 77 1 with the following integral equation 



y(t) = X(i, U) J* PX-^s, U)f (y(s), s)ds - X(t, U) J*' (Id - P)X- 1 (s, U 



){{y(s),s)ds, (80) 



where P is the projection operator associated with the exponential dichotomy (75 I and X is the 



fundamental solution matrix associated with the linear part of (77 1. It can be easily checked 



that the solution of ( 80 1 represents the only solution of ( 77 1 which does not exhibit exponential 



growth or decay within /. Furthermore, using very similar techniques to those employed in [50], 



it can be shown that, for given x(t), the solution of (801 is finite-time hyperbolic and unique 
on the time interval / provided that 



\v(x(t),t) - d x v(x(t),i)(x(i) -x(t))-x(t)\ 



< 00, 



v tei, 



(81) 



and 



\\d*v{x(t),t) -a x v(i(t),t)||oo < I— + -gj , V tel. (82) 
The constants K, L, a, ft are associated with the exponential dichotomy of the linear part of 



(77l (cf. Definition A.3) 
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Definition A. 6 (Distinguished, Finite-time Hyperbolic Trajectory). Let x(t) be a finite-time 
hyperbolic path which does not have an exponential component within I. A trajectory -y(t) 
of the system ( 69 1 is called a Finite-time Distinguished Hyperbolic Trajectory if it can be rep- 
resented as 7(t) = y(t) + x(t) where y(t) satisfies the integral equation (80 1 subject to the 
conditions (811 and (821, and the path x is frozen-time hyperbolic (cf Definition A-5 1 . 

Remarks. Two issues are worth noticing here: 

(i) The frozen-time hyperbolic path x(t) used in Definition | A . 6| can be given, in particular, 
by a path of Instantaneous stagnation points (ISPs) which are frozen-time hyperbolic (cf. Def- 
inition A. 5 I. Given the velocity field v : R n x/^ R n , a path of ISPs is given by a continuous 



curve, Xi 3 p(t), such that 

v(x isp (t),t) = o, tei, (83) 

where T C / is a time interval within which the Jacobian, d x v(Xi 3p (t),t), does not vanish, as 
required by the Implicit Function Theorem for the existence of a solution to (83 1. 

(ii) The notion of a Distinguished, Finite-time hyperbolic trajectory is, in general, non- 
unique on any finite time (or semi-finite) interval. 

Definition A. 7 (Rate- of- Strain tensor). The symmetric part 

|[A 7 (t)+A 7 (t) T ], 



of A~f(t) = d^v(-y(t),t) is called the rate of strain tensor. 



(84) 



The rate of strain tensor describes the growth or decay of solutions £(t) of the linearised 



system (701. This can be seen by directly evaluating d||£(t)|| /dt, i.e. 



(85) 



where (•, •} denotes the canonical inner product on R n , which induces the norm ||£|| = y (£,£)■ 
Thus, if S(t) is negative definite, all solutions of the linearised system are strictly monotonically 
decaying (in the sense of their norm) to the trivial solution. When S(t) is positive definite, 
all solutions of the linearised system are strictly monotonically growing (in the sense of their 
norm). If the strain tensor is indefinite or semi-definite one can define the following set 

Definition A. 8 (Zero-Strain set, cf J&Ef). The set 

(86) 



Zy(t) = {£ & R 2 ■. ($,Sfim}=o}, 



is called the zero-strain set associated with linearisation about 7(i). 

Definition A. 9 (Strain Acceleration Tensor (or Cotter-Rivlin rate of S tensor)). The time 
dependent operator 

d 



M- ( (t) 



dt 



S 7 (t) + S 7 (t)A 7 (t) + A^ty S 7 (i), 



(87) 



is called the strain acceleration set associated with linearisation about "y(t). The restriction of 
My to the zero-strain set is denoted by MS . 

The strain acceleration tensor is associated with the second derivative of i.e. 



dt 2 



dt 



(t(t),s(t)m = (m^(t)m- 



The following extends the dynamic EPH partition of 1R 2 introduced in Haller [38] and 
generalised to a compressible flow setting in Due & Siegmund 28 : 

Definition A. 10 (Dynamic partition of R 2 ). Consider the extended phase space, M 2 x I, 
associated with the flow induced by ( |69[ ). For each t £ I one can define the following sets 



(i) 


Attracting region: 


Ait) 


= {x € R 2 


S x (t) 


is negative definite}, 


(ii) 


Repelling region: 


nit) 


= {xeM 2 


S x (t) 


is positive definite}, 


(iii) 


Elliptic region: 


£(t) 


= {x e R 2 : 


S x (t) 


is indefinite, M x (t) is indefinite}, 


(iv) 


Hyperbolic region: 


H(t)-- 


= {x £ R 2 : 


S x (t) 


is indefinite, M x (t) is positive definite}, 


(iii) 


Quasi-hyperbolic region: 


Qit) = 


--{xeR 2 : 


S x (t) 


s indefinite, M x (t) is negative definite}, 



35 



(iii) Degenerate region: V(t) = R 2 \[A(t) U Tl(t) U £(t) U H(t) U Q(t)}}. 



Definition A. 11 ( Fin ite-time hyperboltcity according to the EPH partition; Holler \38f). As- 
sume that n — 2 in ( 69 1 and that the velocity field satisfies Vw = 0. A trajectory y(t) : I — > R 2 



of (69 1 is called finite-time hyperbolic on the interval / if 

(i) 7(i) intersects T>(I) at isolated points. 

(ii) If Is denotes a time interval that the trajectory spends in £ (I), then 



jf v/2|S 7 (t)|di< |, 



(89) 



where \S\ = ,/E,, =1 |S 



The condition (ii) implies that if *y(t) is finite-time hyperbolic, its finite-time stable and unstable 
manifolds of are non-empty. See [38] for details. See also Appendix IB] and |28l Theorem 42]. 



B On the choice of the initial material segment in nu- 
merical computations of stable and unstable manifolds 
of (finite-time) hyperbolic trajectories. 

We briefly discuss here the problem of approximating stable and unstable manifolds of flow 



trajectories which are finite-time hyperbolic (cf Definition A.4\. 

Consider the linearisation ( |70| l of the dynamical system ( |69| l about a system trajectory (so 
that j(t) — f{~f(t), t) for t £ I). In such a case the the stability properties of the trivial solution 
£(t) = of (70 1 correspond to the linear stability properties of 7(i) in (69 1. 

As already noted in j |3.2.2| if the system ( |69[ ) is only known (or defined) on a bounded 
interval J C JR, it is not possible to define the stable and unstable manifolds of £(t) = in the 
traditional 'inifinite-time' sense even if £(t) =0 is hyperbolic (in the inifinite-time sense) for 
the system ( 70 1 considered on / = R. 

However, if £(t) = is finite-time hyperbolic on 7, one can define (cf [28 ) the following two 
flow-invariant, 'stable' and 'unstable' sets: The finite-time stable set of 7(t) = on / is given 
by 

d 



dm 



|X(m,t)£ t || < 0, V me I 



and the finite-time unstable set of £(t) — on / is defined, for t G /, as 



[ 7 ] = \^ t ,t)eR 2 xl . 



-I 

dm 1 



X(m,t)& || > 0, Vmel 



(90) 



(91) 



where X is the fundamental solution matrix associated with (69 1 and || ■ || is the norm induced 
by the canonical inner product on R 2 , i.e. ||x|| = \J (x,x). The instantaneous geometry of (901 
and (|9T| is given by 



w?t7](*) = U*e B 



and 



*ift](t) = {t t eMt 



^1 

dm 1 

^1 
dm 1 



X(m,i)f t || < 0, Vm G I 



X(m,t)$ t || > 0, VmG/ 



(92) 



(93) 



referred to as i-fibres. In contrast to the classical (time asymptotic) definition of stable and un- 
stable manifolds, the finite-time counterparts, W™ and Wf , have the dimension of the extended 
phase space (rather than a lower dimension) and their t-fibres are open sets in R n . In such 
a case, a common approach used in the invariant-manifold Lagrangian transport analysis is to 
choose (non-unique) segments of initial conditions of length a< 1, iXf a and & f a , containing the 
trivial solution of the linearised system and follow their forwards and backward time evolution. 



We show below (cf Proposition B.2 l how to choose the (non-unique) material segments in such 



a way that they are contained in, respectively, the finite-time stable and unstable manifolds. 



3G 



Recall first that the trivial solution £(t) = of the linearised equation corresponds to 7(i) 
of (691. Thus, if 7(t) G T~L(t) then also G 7i(t). Furthermore, if £ is fi nite-tim e hyperbolic 
on 7, the (symmetric) rate of strain tensor, S(t), is indefinite (cf Definition A. 10 1 on I so that, 
for each t £ I, the zero-strain set contains two orthogonal lines and is given by 

Zf(t) = {zi,Z2 £ R 2 '■ {zi(t),za(t))=0, (*i(t),S 7 (t)z 2 (*)) =0}. (94) 
We now define a subset of nondecreasing solutions at t as 



* + (t) = U* G JR 



dm 



|X(m,t)f t 



>0k 



a subset of nonincreasing solutions at t in as 



*"(*) = {it e R 



At 



|X(m,t)& 



II < 



(95) 



(96) 



so that * + (t) n * (i) = Z~ t {t)- Moreover, for £ H(t) the restriction of the strain 

acceleration tensor to Z~/(t), M 7 (t), is positive definite, i.e. (£i(t), M(t)£i(i)) > and 
(fa(<),JWr(t){2(t)) >0which, based on ( 88 1, implies that solutions, £(t, ft* , t*), ft* G Z~/(t*), of 



(70 1 cross the zero strain set Z 1 {t*) at t = t* from the region of decreasing norm to the region 
of increasing norm. 



Proposition B.l. Consider a trajectory j(t) of (691 and the corresponding trivial solution of 
the linearised system | |70[ ) on I = [t a ,tfi] with the fundamental solution matrix X(t, t a ). The 
finite-time unstable set, Wf py(i)] and the finite-time stable set, Wf \y(t)] are invariant under the 
action of X(t,t ). Moreover, if 7(f) G H(t) fort G I, the set = {£ G I? 2 : 3 t G J, £ G * + (t)} 
is forward-time invariant and the set tyj — {£ G M 2 : 3 t G /, £ G $?~(t)} is backward time 
invariant. In particular, Wf [7]^) = * + (ta) and Wf ['y](t b ) = *~(*b). 

Proof. The invariance of Wf [7(f)] and Wf [7(f)], as well as the forward-time invariance of ^ + (t) 
and the backward-time invariance of ^~(t), was discussed in [281 cf Remark 23, Theorem 44]. 
In order to show that Wf [7]^) = ^ ,+ {ta) we appeal to the forward invariance of ^ + {t) under 
the action of X(i, t a ). 

Assume first that the opposite holds, i.e. that ft* G * + (t*) andthat£(i**,£ t »,i*) $ * + (t") 
for t* < £**, t*,t** G 7. Due to continuity of f(i), the trajectory has to cross the zero strain 
set at some time t* < t x < t** which requires that £(i x ,ft*,£*) G Z 1 (t x ) and 



At 2 



l|f(*,ft*,0 



<f(t x ),M(i x )£(t x ))<0, f t * G *+(£*), 



(97) 



which contradicts the fact that if -y(t) G Ti.(t) for t £ I, M(t) is positive definite on Z 1 {t) for 
f G J. Consequently, if j(t) £ H(t) and f t « G *+(**), then f(i) G * + (i) for t > t*, t,f G /, 
which implies that ty + is forward-time invariant on I. Note that "if + is not backward time 
invariant. In order to see this, it is sufficient to consider trajectories crossing the zero strain 
set, Z y (t*), at t* £ (t a ,tb\. Since 9* + (t*) = Z 1 (t*), any trajectory f(t,f t ., i*),f t « G Z 1 (t*) 
leaves for t < t* in backward time. We finally not that the set ty + (t a ) is invariant under 
the action of X(t,t a ), which implies that ty + (t a ) C Wf[7](t ). However, based on definitions 
(95} and Q is is clear that Wf[y](t) C * + (t) which implies that * + (i a ) = Wf [y}(t a ). 

Similar procedure can be used in backward time to show that is backward-time invariant 



on I. 



□ 



Proposition B.2. Consider the linearised flow (70 1 over the time interval I so that the trivial 
solution is finite-time hyperbolic on I (in the sense that j(t) £ Ti for t £ I). If the material 



U? a = {a; G R 2 : x = (3S + (t a 



and 



/3e 



&t b = {X £ R 2 : X = PS~(t b ) 13 £ 



a a. 
~2' 2 

a a 
"2' 2 1 



C R}, 
C R}, 



(98) 



(99) 



where S + (t) and S (t) are the eigenvectors of the rate of strain tensor, S(t), corresponding to 
the eigenvalues s+(i) > 0, s~(t) < 0. Then, U? a C W? a {f = 0} and 6? a C Wfjf = 0}. 
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Proof. For any point ut a = f3S + (t a ) G ii" a , \(3\ ^ a/2, we have 

(u ta ,S(t a )ut a ) > 0, (100) 

which implies that u ta G ^ + (t a ) = W/py](i a ). The invariance of W"[7] implies that 
t(t,U ta ,t a ) G W?py](i) for t G J. Similarly, for any point St„ = /3S (t a ) G ©£, j/3j sC a/2, we 
have 

K.sfttK) <o, (ioi) 

which implies that St b G *~(tb) = Wfpy](t{,). The invariance of Wfpy] implies that 

e(*,*t»,*»)ewy[7](*) forte/. □ 

Note finally that, due to the the embedding property of finite-time stable and unstable 
manifolds (see |28l Theorem 37, p. 659]), the stable and unstable manifolds of £(t), for two 
time intervals /, J, such that I C J, satisfy the following 

W! C Wj and W? (ZWj. (102) 

Thus, the effect of the non-unique choice of the initial material segments diminishes with the 
length of the considered time interval, provided that the considered trajectory is finite-time 
hyperbolic on I. 
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